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ABSTRACT 


An investigation is ma<3e of the problem of predicting 
the attitude of satellites under the influence of external 
disturbing torques. The attitude dynamics are first 
expressed in a perturbation f ojrmulation which is then 
solved by the multiple scales approach. The independent 
variable, time, is extended into new scales, fast, slow, 
etc., and the integration is carried out separately in 

I 

the new variables. The rapid and slow aspects of the 
dynamics are thus systematically separated, resulting 
in a more rapid computer implementation. The theory is 
applied to two different satellite configurations, rigid 
body and dual spin, each of which may have an asymmetric 
mass distribution. The disturbing torques considered 
are gravity gradient and geomagnetic. A comparison 
with conventional numerical integration shows that our 
approach is faster by an order of magnitude. 

Finally, as multiple time scales approach separates 
slow and fast behaviors of satellite attitude motion, 
this property is used for the design of an attitude 
control device. A nutation damping control loop, using 
the geomagnetic torque for an earth pointing dual spin 
satellite, is designed in terms of the slow equation. 
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CHAPTER 1 


INTRODUCTION 


1.1 General Background 

The problem of predicting a satellite attitude motion, 
under the influence of its environmental torques , is of 
fundaraental importance to many problems in space research. 
An example is the determination of required control torque 
as well as the amount of fuel or energy for the satellite 
attitude control devices. Similarly, a better prediction 
of the satellite attitude motion can be helpful in yielding 
more accurate data for many onboard experiments, such as 
the measurement of the geomagnetic field or the upper 
atmosphere density etc., which depend on the satellite 
orientation. 

Yet, the problem of satellite attitude prediction 
is still one of the more difficult problems confronting 
space engineers today. Mathematically, the problem con- 
sists of integrating a set of non-linear differential 
equations with given initial conditions , such that the 
satellite attitude motion can be found as functions of 
time. However, the process of integrating these 
equations by a direct numerical method for long time inter- 
vals, such as hours, days (which could be even months or 
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years ) , is practically prohibited for reasons of compu- 
tational cost and possible propagation of numerical 
round-off and truncation errors. On the other hand, it 
is even more difficult, if it is possible at all, to 
have an exact analytic solution of the problem, because 
of the non-linearity and the existence of various external 
disturbing torques in each circumstance. 

A reasonable alternative approach to the above problem 
seems to be to apply an asymptotic technique for yielding an 
approximate solution. The purpose of this approach is 
to reduce the computational effort in the task of attitude 
prediction for long intervals, at the cost of introducing 
some asymptotic approximation errors. Meanwhile, the 
asymptotic approximate solution has to be numerically 
implemented in order to make it capable of handling a 
broad class of situations. 

We found the problem is interesting and challenging 
in two ways. First , because it is basic, the results 
may have many applications. Second, the problem is very 
complicated, even an approximate approach is difficult 
both from analytic and numerical points of view. 


'OMGHAU'PAUBl^ 
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1.2 Problem Description 

The subject of this thesis is the prediction of 
satellite attitude motion under the influence of various 
disturbing torques. The main objective is to formulate 
a fast and accurate way of simulating the attitude rota- 
tional dynamics in terms of the angular velocity and 
Euler parameters as functions of time. The formulation 
has to be general, it must be able to handle any orbit, 
initial conditions, or satellite mass distribution. Fur- 
thermore, it must predict the long term secular effects 
and / or the complete attitude rotational motion, depending 
on the requirement. Because of this built-in generality 
it is intended that the program can be used as a design 
tool for many practical .space engineering designs. To 
achieve this desired end the problem is first expressed 
as an Encke formulation. Then, the multiple time scales 
(MTS) technique is applied to obtain a uniformly valid 
asymptotic approximate solution to first order for the 
perturbed attitude dyneimics , 

Two different satellite configurations are considered, 
a rigid body satellite and a dual spin satellite, each of 
which may have an asymmetric mass distribution. In the 
latter case, it is assumed that the satellite contains a 
single fly wheel, mounted along one of the satellite 
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body-principal-axes, to stabilize the satellite attitude 
motion. These models are considered typical of many 
classes of satellites in operation today. The disturbing 
torques considered in this dissertation are the gravity 
gradient and the geomagnetic torques . For a high-orbit 
earth satellite these two torques are at least a hundred 
times bigger than any other possible disturbance, though, 
of course, there would be no difficulty inserting models 
of other perturbations. 

Both the gravity gradient and the geomagnetic torques 
depend on the position as well as the attitude of the sate- 
llite with respect to the earth. Therefore, the orbital 
and attitude motion are slowly mixed by the actions of 
these disturbances. Hov/ever, the attitude motion of the 
vehicle about its center of mass could occur at a much 
faster rate than the motion of the vehicle in orbit around 
the earth. Directly integrating this mixed motion, fast 
and slow together, is very inefficient in terms of com- 
puter time. However, realizing that there are these 
different rates, then the ratio of the averaged orbital 
angular velocity, to the averaged attitude angular velo- 
city or equivalently the ratio of the orbital and attitude 
frequencies (a small parameter denoted £ ) may be used in 
the MTS technique to separate the slow orbital motion 
from the fast attitude motion. In this way the original 

11 
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dynamics are replaced by two differential equations in 
terms of a slow and a fast time scale respectively. 

The secular effect as i^ell as the orbit-attitude coupling is 
then given by the equation in the slow time scale, while- 
the non-biased oscillatory motion is given by the second 
equation in terms of the fast time scale. In addition, 
a method for handling the resonance problem is also discussed. 

In some situations the slow equation for the secular 
effects can be useful in the design of an attitude control 
system. The vehicle environment torques, if properly 
used, can be harnessed as a control force. However, to 
design such a control system it is often found that the 
control force is much too small to analyze the problem 
in the usual way. In fact, the design is facilitated in 
terms of the equation of the slow variable because only 
the long term secular motions can be affected. This 
application is demonstrated by mean of a nutation damping 
control loop using the geomagnetic torque. 
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1,3 His-torical Review and Literature Survey 

The problem of attitude dynamics of an artificial 
satellite, rigid body as well as dual spin case, is 
closely related to the branch of mechanics of rigid body 
rotational motion. The subject is regarded as one of 
the oldest branches of science, starting from middle of 
the eighteenth century. Since then it has interested 
many brilliant minds for generations. The literature 
in this area, therefore, is rich and vast. Thus, we have 
to focus our attention on only those areas which are 
immediately related to this research. 

The classical approach to the rotational dynamics 
mainly seeks the analytic solutions and their geometric 
interpretations. By this approach, many important and 
elegant results have been obtained. Among them, the 
'Poinsot construction* of L, Poinsot (11] , gives a geo- 
metrical representation of the rigid body rotational 
motion. The Suler-Poinsot problem, for a torque-free 
motion, was first solved by G, Kirchhof f (12] by means of 
Jacobian elliptic functions. On the other hand, F, 

Klein and A. Sommerfeld [13] , formulated the same problem 
in terms of the singularity-free Euler symmetric para- 
meters and gave the solution. Recently, H. Morton, J, 
Junkins and others [15] solved the equations of Euler 
symmetric parameters again by introducing a set of complex 
orientation parameters, Kirchhoff's solution as well as 
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the solution of Euler symmetric parameters by H. Morton etc. ‘ 
play an important role- as a reference trajectory- in our 
study. 

This approach, however, because of its nature can 
not handle a general problem for various situations. This 
difficulty is substantial for the case of an artificial 
satellite. A more flexible alternative, widely applied in 
the engineering world, is the asymptotic technique for 
evaluating an approximate solution. Among them, the 
averaging method by N. N. Bogoliubov and Y. Mitropolsky 
[14J is the most commonly used. For example Holland and 
Sperling [16] have used the averaging method for estimating 
the slow variational motion of the satellite angular 
momentiam vector under the influence of gravity gradient 
torque and Beletskii [17] formulated perturbation 
equations using the osculating elements for a dynamically 
symmetric satellite, P, L. Chernous'ko [18] derived 
the equations of variation of parameters for angular 
momentum vector and the rotational kinetic energy for 
an asymmetric satellite. However, the averaging method 
is most easily applied for a problem which normally has 
a set of constant parameters, such that the slow variational 
behavior of these parameters can be established in a 
perturbed situation. For example in a simple harmonic 
oscillator, the frequency and amplitude are two para- 
meters which characterize the dynamics described by a 


14 



second order differential equation. Unfortunately, 
rotational motion in general, does not immediately lead 
to a complete set of similar parameters. Although, it 
has constant angular momentum vector and kinetic energy 
as parameters it is a six-dimensional problem. Besides, 
an elliptic integral is involed in its solution. Never- 
theless, this difficulty can be overcome by casting the 
problem in a Hamilton- Jacobi form, from which a variation- 
of-parameter formulation can be derived in terms of 
Jacobi elements. This approach is reflected in the works 
of Hitzl and Breakwell [19] Cochran[20], Pringle[21] etc. 

Our dissertation is different from the others 
mainly in three aspects. First, it is a new approach, 
using the multiple time-scales method [1-7] with the 
Encke perturbation formulation [22] , for predicting the 
complete satellite attitude motion without involving the 
Hamilton- Jacobi equation. Second, we are interested in 
the secular effect of the disturbing torques as well 
as the non-secular oscillatory effect. By combining 
them, we have the complete solution. Further we know 
that the former gives the long term behavior and the 
latter indicates the high-frequency motion of the sate- 
llite attitude dynamics. Third, our immediate objective 
is numerically oriented for saving computer time. Thus 
the difficulties we encounter could be analytical as 
well as numerical . 
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1.4 Arrangement of the Dissertation 

Chapter 2 reviews the multiple time scales asymptotic 
technique - a basic tool in this research. Two examples 
are used for illustrating the fundamental procedure, one 
represents the secular type of almost-linear problem 
and the other represents the singular type of slowly 
time-varying linear system. 

Chapter 3 develops the asymptotic solution to the 
attitude motion of a rigid body satellite under the influ- 
ence of known small external torques. It shows that, 
the original equations of attitude dynamics can be repre- 
sented by two separate equations - one describing the 
slow secular effects, and the other describing the fast 
oscillatory motion. The latter can be analytically 
evaluated. Numerical simulation using this approach is 
also presented for the class of rigid body satellites 
under the influence of gravity gradient and geomagnetic 
torques . 

In chapter 4, the previous results are extended 
to the case of dual spin satellite, in which .a fly-wheel 
is mounted onboard. Two sets of numerical simulations, 
one for a dual-spin satellite in the earth gravity 
gradient field and the other influenced by the geomagnetic 
field, are given. 
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Chapter 5 represents an application of the fact that 
MTS method separates the slow and fast behaviors of a satel- 
lite attitude motion. We demonstrate that the slow equation, 
which describes the secular effects , can be useful in the 
design of a satellite attitude control system, A 
nutation damping feedback control loop, using the geomag- 
netic torque for an earth pointing dual-spin satellite, 
is designed in terms of the slow equation. 

In chapter 6 the conclusions drawn from the results 
of this study are summarized, and some suggestions for 
future research are listed. 




CHAPTER 2 


REVIEW OF MULTIPLE TIME SCALES (MTS) METHOD 
2.1 Introduction 

multiple time scales (MTS) method is one of the 
relatively newly developed asymptotic techniques. It 
enables us to develop approximate solutions to some 
complicated problems involving a small parameter £ , 

when the exact solutions are difficult, if not impossible, 
to find. The basic concept of MTS method is to extend 
the independent variable, usually time, into multi-dimen- 
sions. They are then used together with the expansion of 
the solution (dependent variable) such that an extra degree 
of freedom is created and the artificial secular terms can 
be removed. Thus a uniformly valid approximate solution 
is obtained {1-7] . 

An unique feature of the MTS method is that it can 
handle secular type as well as singular type of perturbation 
problems in a unified approach. By this method, the fast 
and slow behaviors of the dynamics are systematically iden- 
tified and separated. The rapid motion is given in terms 
of a fast time scale and the slow motion in terms of a 
slow' time scale, each of which, in most cases, has a 
meaningful physical explanation. A comprehensive refer- 
ence on this subject is by Ramnath[3]. The textbook by 
Nayfeh [7] has also been found informative. 
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2.2 MTS and Secular Type of Problems 


A secular perturbation problem is one in which 
the nonuniformity in a direct expansion occurs for large 
values of the independent variable. We consider systems 
with a small forcing term. The forcing term changes 
the dynamics gradually and has no appreciable effect in 
a short time. However, the long time secular effect of 
the small forcing term may significantly influence the 
overall behavior of the dynamics. From a mathematical 
point of view, a secular type of problem has a singularity 
at infinty in the time domain. 

Since perturbation problems and the asymptotic 
technique for solving them can be most easily understood 
by solving a demonstration case, let us consider a simple 
example of a slowly damped linear oscillator [7] , 

X + X = - Z £ X y 0<£<-<! {2.2.1) 

where £ is a small parameter. The simplicity of the 
forcing term ~2JtX allows us to interpret the approximate 
solution, developed later. The exact solution is available, 
but the generality of our asymptotic approach will not be 
lost in spite of the simple form of the forcing term. 

We first solve the problem by Poincare type of . 

direct expansion method [33] , such that difficulties of 

non-uniformity and secular terms can be illustrated. 

Then the same problem is solved by MTS method, which 
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yields a uniformly valid asymptotic solution to first 
order . 

We expand X into an asymptotic series in £ t 
0C{.±) ~ Sc d- ‘ ' 

= 4 if- ^)C2,i±) ,9 0 ^ 


An asymptotic series is defined [9] as one in which 
the magnitude of each term is at least one order less than 
its previous one I I / 1 £>^-1 1 = I £*^Xn ( /| j 

Therefore Sx. decreases rapidly as the index i increases. 
This simple fact allows us to approximate the solution by 
calculating only a few leading terms in the series expansion. 

Substituting (2.1.2) into (2.1.1) and equating the 
coefficients of like powers of £ , we have: • 


Xo + ^0=0 



X, -+ X , = -:£-Xo 


(2.2.3) 


(2.2.4) 


-cXc 
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The solution for in equation (2.2.3) is 

Tdj, = iX Cc-a-^i: + b ^ (2.2.5) 

Where 'a' and ’b’ are two constants. Substituting 
into eq, (2,2.4) of Xj 

+X, = (-0^ ^ (2.2.6) 

lui-kK X. C. Xi( 0 )= X, {D) = o 

The solution is I 

X I s:. — o^'k Ccx ^ 1o ^ Awn, ^ -4- A A. 

' (2.2.7) 

The approximation of x(t) up to first order is 
therefore: 

XoC*)4 ‘LX/tt) 

— ( (X Ccsi± 4 . b-4^4) 

•+1 (X. A “ \pi: At><^A -f- 

( 2 . 2 . 8 ) 

Above approximation, however, is a poor one because 
of the occurrence of two terms — Cea A and — b/t 

which approach infinity a.s . They are referred 

to as secular terms. We know the true x(t) has to be 
bounded and asymptotically decaying, for x(t) is a . 
damped' harmonic, oscillator. The secular terms 
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make the series expansion 

eq. (2,2.8) a non-asymptotic one, since »■ co 

3o 

■as ^ -^oo . In the process of finding a solution by series 
expansion, there is no guarantee that the higher order terms 
can be ignored in a non-asymptotic series expansion. On the 
other hand, if an asymptotic series is truncated, the error 
due to the ignored higher order terms will be uniformly 
bounded in a sense that j error j / (approximate solution) 

. An approximate solution is said to be 'uniformly 
valid' if its error is uniformly bounded in the interval of 
interest. We see that the loss of accuracy 
by straightforward Poincare type of expansion is due to the 
occurrence of the secular terms and therefore the approxi- 
mation is not uniformly valid. 

In the following, the same problem will be studied in 
the context of multiple time scale method, which yields a 
uniformly valid solution. 

For convenience, we rewrite the dynamics 




( 2 , 2 . 9 ) 


The solution is first expanded into an asymptotic 

series of t , same as before 


Xti) := Xo(*>-+ £ -t 


( 2 . 2 . 10 ) 


The concept of extension is then invoked, which 
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expands the domain of the independent variable into a space 
of many independent variables. These new independent 
variables as well as the new terms that arise due to 
extension are then so chosen that the non-uniformities 
of direct perturbation can be eliminated [3] . 

Let 







( 2 . 2 . 11 ) 


For an almost linear problem, as the one we are studying 

% 

jiXc 


The new time scales, in general, can be non-linear 
as well as complex, v;hich depend upon the nature of the 
problem [4] . However, for this particular problem, a set of 
simple linear time scales works just as well. 

The time derivative can be extended in the 
space by partial derivatives as follows. 


jL 

di: 


'^'^2. dd: 






+ 1 


a-^i 










23 



And 


JL 




+ XI 








^ V 3Xt3Xj ai=,^ ^ 


Substituting equations (2.2.12), (2.2.13) and (2.2.10) 
into equation (2.2.9), and equating coefficients of like 
powers of £- , we have 


: 


~b 2Cq 

Tf- 


4 '>ic = O 


(2.2.14) 


a T/ ' 


■=-2 


9 X 


2 




; (2.2.15) 
_cXc. 

The original equation has been replaced by a set of 
partial differential equations. The solution for ~>Lo(^) 
from (2.2.14) is 


Xo — (X, , Xi - ) (aZo) -h h (X/ , Ci ■• ) (-A Zo) 


(2.2.16) 


where a,b are functions of ,...etc., and are 

yet to be determined. Substitute from (2.2.16) into 
(2.2.15), and solve for . 


<> X j 


4 %j 


4 jtoCp(-viTo> ( 4 bA 44 X 


(2.2.17) 
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Equation (2.2,17) represents a harmonic oscillator 
driven by an external sinusoidal function. The 

could become unlimited^ since the external force 
has the same natural frequency as the system itself. 

In order to have \Xf/Xo\ bounded, the terms to the right 
hand side of the equal sign in (2.2.17) have to be 

set to zero, which will result in a bounded Xi(^) . Note 
that, this is possible, because there is a freedom of 
selecting functions 'a' and 'b* from the extension of the 
independent variable. By doing so , we have. 



“P 2, ' 

- 0 


Z 'C^ 



+ Z. 


(2.2.18) 


ov 


f (X ~ do (~ h ) 

(2.2.19) 

b = 

where do and bg are two constants. Combining 
(2.2.16) and (2. 2. 19), the approximate solution for X(.±) 
up to first order of t by MTS method is 

X(t) = Xo + 

= bo • 

( 2 . 2 . 20 ) 
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The exact solution for (2,2,9) can be obtained, 
which is 

4- bo ^pC-^rlT - (2.2.21) 

The error/solution ratio in this case is 


{ er\ror 


X exftjtt ~ X~cipprcj( 

Solution 


XjLyac-t 


0 


It is interesting to note that by the MTS method, 
we have replaced the original dynamics (2.2.9) by an 
equation (2.2,14) in the fast time scale ^ 

'h^^c 


5 


+ Xc ••= o 


(2.2.14) 


and two slow equations in the slow time scale . 


ax, 


->r - Q 


b b 


I 


+ b = o 


(2.2.18) 


The fast equation gives the undisturbed oscillatory 
motion and the slow equations represent the slow varia- 
tional change of the amplitude of the oscillation caused 
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by the damping term 


2.3 MTS Method and Singular Perturbation Problems 

• There is a class of perturbation problems in which 
the behavior of the reduced system - by setting £ 
equal to zero - could be dramatically different from the 
original. This phenomenon occurs because the reduced 
system, described by a lower order differential equation 
can not satisfy the given boundary conditions in general. 
We call this kind of problem a singular perturbation 
problem. 

Singular perturbation problems have played an 
important role in the engineering field, most notably, 
in fluid mechanics, for eg. the boundary layer theory. 

This problem v/as solved by introducing the inner (Prandtl' 
boundary layer) and outer expansion [32] . However, 
the same problem also can be solved in a more straight- 
forward approach by the MTS method. This approach was 
first noted in the paper by Ramnath [4] in studying the 
behavior of a slowly time-variant linear system by 
employing non-linear time scales. In the following, let 
us use an example, which is adopted from [3] , for demon- 
stration. 
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Consider a second order singular perturbation problem ' 


y + uj itj) y = o 


(2.3.1) 


Where o< £« 1 , is a constant small parameter. 

Expand the time domain into two -dime ns ions 


X ~ ' I-} 3 

and define Tc^ 'tj as -follows 




rr = J- 




(2.3.2) 


where is yet to be determined. 

j ^ 

The time derivative - 7 ^ and — rrr 

dt dx*’ 


can be extended as 


<K-k ^ t dx, ‘ 5x^ 

I ^ I 3^ ^ (2.3.3) 

4 _! ^ ^ ^ 3^ 

£, 3X; 3x/ 

Substituting (2.3.3) into (2.3.1) and separating terms 
according to the pover of S , we will have the set of 
equations , 


X- 
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£‘*-- 


:ir 

+ i\Ji-Ce) y = o 

(2.3.4) 



(2.3.5) 


= 0 

(2.3.6) 

3 V 




By assuming that ^ solution in the form 


^ "r,) = ‘sCCTo ) -a/x/> (x, ) 

substitution into (2.3,4), yields 


-+ tJC'Tc) = 


or 


= X OJ^ 

= "A aJ ^ 


(2.3.7) 


(2.3.8) 


Similarly put (2.3.7) into (2.3.5), we have 

- i -A.. c( -h Z ~ O (2.3.9) 

3 To ^'^0 

oC - ^ ^ (2.3.10) 

The approximate solution up to first order can be 
constructed as 
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-/ -I 

y(« = c; iir J 

zi ' J. — ' - 1 

~~ UJ^ ['^ -{r -C^y3 [ <^/(J 

(2.3.11) 

We obtain (2.3.11) from (2.3.8) and (2,3.10). 

Note that in our approximation the frequency variation 
is described on the ■~r, scale and the amplitude variation 
on the scale. The success of this approach depends 
on the proper choice of the nonlinear clock. While in 
the past this choice was made on intuitive grounds, 
recent work [ 2. ] has been directed towards a systematic 
determation of the clocks . Ramnath [ 3 ] has shown that 
the best nonlinear clocks can be determined purely in 
a deductive manner, from the governing equations of the 
system. With a judicious choice of scales the accuracy 
of the asymptotic approximation is assured, A detailed' 
error analysis of the approximation was given by Ramnath 
[ 3 ] . These questions are beyond the scope of the 
present effort and reference [3 ] may be consulted for 

c 

more information. 
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CHAPTER 3 


PREDICTION OF ATTITUDE MOTION FOR A RIGID BODY SATELLITE 


Introduction 

In this chapter, a multiple time scales asymptotic 
technique is applied for the prediction of a rigid body 
satellite attitude motion disturbed by a small external 
torque. 

The attitude dynamics of a satellite, described in 
terms of the Euler's equations and Euler symmetric para- 
meters, are first perturbed into an Encke formulation, 
in which the torque-free case is considered as a nominal 
solution. The multiple time scales technique is then 
used for the separation of the fast attitude motion from 
the slow orbital motion in an approximate, but asymptotic 
way. Thereby, the original dynamics can be replaced by 
two sets of partial differential equations in terms of a 
slow and a fast time scale. The long-term secular effects 
due to the disturbing torque are given by the equations 
in the slow time scale, which operate at the same rate as 
the orbital motion. A non-biased oscillatory motion is 
given by the second set of equations in terras of the fast 
time scale, which basically describes the vehicle attitude 
oscillatory motion. These fast and slow motions, combined. 
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give us a first order asymptotic solution to the Encke’s 
perturbational equation, whichj therefore j can be regarded 
as the second order asymptotic solution to the original 
satellite attitude dynamics. 

Finally, the fast non-biased oscillatory motion 
can be analytically evaluated if the external torques 
are not explicitly functions of time. Thus numerical 
simulation of a satellite rotational motion by this 
new approach requires only the integration of the slow 
equation which can be done with a large integration 
time step. This fact leads to a significant saving of 
computer time as compared to a direct numerical inte- 
gration, Two examples, one with gravity gradient torque 
the other with geomagnetic torque, are demonstrated in 
section. 3.6 and 3.7. 


3.2 Rigid Body Rotational Dynamics 
(A) . Euler's Equations 

Newton's second law for rigid body rotational motion 
in an inertial frame can be written as 






(3.2.1) 


Where H and N are the angular momentum and the ' 
external torque. By Coriolis law, the motion can be 
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expressed in any moving frame *b' , as: 

ol H !'“i b ^ ^ ^ in b 

4- Kb {3.2.2 

(A X 

Where is the angular velocity of the 'b' 

frame with respect to the inertial frame. In case the 
'b' frame is selected to coincide with the body fixed 
principal axes (x,y,z), then the angular momentum can be 
written as 



(3.2.3) 

where f ly $ are moments of inertia about x, y, z axes. 
Combining (3.2.2) and (3.2.3) we have Euler's equations 

Ix - ( Y 

I^) ( 3 . 2 . 4 ) 

In vector notation, they are ; • 
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(3.2.5) 


VO -+ ( IV ;< ) uj = M 

Euler's equations give the angular velocity of a 
rigid body with respect to an inertial space though this 
angular velocity is expressed in the instantaneous body 
fixed principal axes [35] . 

(B) . Euler Symiaetric Parameters 

The role of Euler syimnetric parameters are 
similar to Euler angles, which define the relative 
orientation between two coordinates. From either of 
them, a transformation matrix can be calculated and a 
vector can be transformed from one coordinate to another 
by pre-multiplying with the transformation matrix. 

However, from an application point of view, there are 
notable differences between Euler symmetric parameters 
and Euler angles. The important ones are listed as follows. 

1. Euler angles (0^ X=:l;'^y3) have order of three j 

* 

whereas Euler symmetric parameters ( ^ Z^o, t -2,3) 

have order of four with one constraint. 

2, are free from singularity, where QZ are 

not. Since $i ^ (5^ ^9^ are the z-x-z rotations, 
in case that 9^~o , one can not distinguish. 
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from 63 , 

3, |5^ propagate by a linear homogenous differential 
equation; $l are by a non-linear differential 
equation. 

4. ex have a clear physical interpretation, -i.e. 
precession, nutation and rotation. can not 
be immediately visualized. 

By considering above differences, we feel that 
Euler symmetric parameters are more suitable for numerical 
computation because they are propagated by a linear 
equation and free from singularity^ even though they add 
one more dimension to the problem. On the other hand, 

Euler angles are easier to understand. In this chapter 
we select Euler symmetric parameters for the satellite 
attitude prediction problem. And in chapter 5, for a 
satellite attitude control system design, we use Euler 
angles . 

The concept of Euler symmetric parameters is based 

' . J 

upon * Euler Theorem', which says that a completely general 

angular displacement of a rigid body can be accomplished 

/\ 

by a single rotation^ about a unit vector 2 / 

A 

where 2 is fixed to both body and reference frames. 

The /5X are then defined as ; 
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z. 


A = 


1 . 2-/ 3 


(3.2.6) 


/3c ^ 

' ~ 2 . 






with a constraint 


3 


X 



(3.2.7) 


a 

If Cyjj is the transformation matrix from the 
reference frame *r* to the frame ’b', 

* 

Then C^j, can be calculated in term of as [35] 





(3.2.8) 


Also satisfy a homogeneous linear differential 
equation [35] ) 



" N 


(/ > 



A 


_ 1 

C> ~ 



(it 


”” 2 






. 





(3.2,9) 
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in X direction, etc 


— b 

where is the component of 

3.3 Euler-Poinsot Problem 

Closed form solutions for a rigid body rotational 
motion with external torque are usually not possible 
except for a few special cases. A particular one, named 
after Euler and Poinsot, is the zero external torque 
case. This is useful here, since the disturbing 
torques acting on a satellite are small, the Euler-Poinsot 
case can be taken as a nominal trajectory. 

It is Kirchhoff [12] who first derived the complete 
analytic solution for Euler's equation ( dV ) in terms of 
time, in which an elliptic integral of the first kind is 
involved. In the following, we will .review the Kirchhoff 's 
solution along with the solution of Euler symmetric 
parameters (6 ) by Morton and Junkins [15] , Also, 

by defining a polhode frequency, we find that the solution 
for j8i can be further simplified, such that it contains 
only periodic functions, 

(A) Kirchhoff *s Solution 

Without external torque, Euler's equations are 





I O 

(3.3.1) 

"3 

uJ 

-Ci|- ix) 


.{3. 3. 2) 


37 



^3 - T;) ) S - 


( 3 . 3 . 3 ) 


Above equations are equivalent to a third order 
homogeneous ordinary differential equation and its solu- 
tion involves three integration constants. 

Multiplying above three equations by ujj( , Wy and 
respectively, and integrating the sum, we obtain one of 
the integration constants for the problem, called T, which 
is the rotational kinetic energy of the system, that is ^ 




( 3 . 3 . 4 ) 


Similarly, by multiplying the three Euler's equations 
by Ix**V • respectively and inte- 

grating the sum, we have H, another integration constant, 
that is the angular momentum of the system. 


+ J/ 


( 3 . 3 . 5 ) 


Having rotational energy T and angular momentum H 
given, a new variable can be defined in terms of time 
by an elliptic integral of the first kind ^ 


xa-Tj = 


c /V / — 


( 3 . 3 . 6 ) 
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where A , x and k are constants, and k is the 
modulus of the elliptic integral. 

Kirchhoff's solution can be written as follows ,* 



a 

( 1' ^ 

(3.3.7) 

ujy =: 

b 

w 

(3.3.8) 

UJ^ = 

C 

C<M ^ 

(3. '3. 9) 


where constants a,b,c,k, and X are defined as follows, 

/ 

(3.3.10a) 








6"- = 


X^ = 




_ ^^xT- H 


^ JxT - 


ly 


Xv ~ - 1 3 


^x- 


:tJxT- H 
H"-2X^ T 


(3.3.10b) 


(3.3.10c) 


(3.3.10d) 


(3.3.10e) 


X = 


^ Jo 


(3.3.10f) 


39 



^0 = * [' 




Signs of a,b,c and )\ should be picked, such that 
they satisfy the equation 


Xy- Iv 


(X h 


The validity of above solution can be proved by 
direct substitution. 

Kirchhoff's solution is less popular than 
Poinsot construction in the engineering world. One 
reason is that it involves an elliptic integral, 
the solution becomes rather unfamiliar to the engineering 
analyst. However, for our long term satellite attitude 
prediction problem, Kirchhoff's solution seems to be a 
powerful tool. 

(B) Solution For Euler Symmetric Parameters 
Euler symmetric parameters satisfy a linear differ- 
ential equation, which relates and [35] 

( 0 -Myt 'W,]( P„) 


di P. . I 

TT- “ 2- 


0 





0 


-OJ.^ 



o 

UJy 
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o 



The constraint is 


A set of complex numbers L^l) can be introduced 

as follows 




•4" 

A. 

od, = 


•+ 

A 

0^2. =‘ 


— 

^ /3i 


^3 

4 

: /3. 


where 

^ =rFT 

and satisfy a constraint of 

oCo oC^ - = j 

In matrix and vector notations , 

are 

_ ^ ~ 

od = A iS 


(3.3.12) and 


(3.3.14) 
(3.3.13) 

(3.3.15) 

(3.3.16) 
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where 


and 





Q _ W/j^ 

^ (Ay 

-“'j 

[^J = 

o 

*^3 




0 



/ 

-UJ^ 

o 

J 


" j O A 

o 


■K 

OX. 0 

- 1 


h - 

1 0 -A 

o 



0 A o 

1 



Substituting (3.3.16) into (3.3.15), we obtain: 


cL 



* 


[tu] 


^-{ - 
/\ c(. 


(3.3.17) 


or 


" 1 
z 

( i 

4 ( 1 i^x - 

] 0^0 

i, = i 

[ A 

Wj, + (i K)^- u;^)^] 

1 

^2 = i 

r-^ 

4 (X ,J^4^5)-^] 

1 



";-4 (,4^y4W3) -^] 

1 ^^3 
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For a torque - free rigid body rotation, the 
angular momentum vector H remains constant. If we assume 
one of our inertial axes (say A^.) pointing in H direction 
and name this particular inertial frame by ’n’, then 


H 


n 


O 

H 

o 


Since 




(3.3.19) 


(3.3.20) 


combining (3.3.19), (3.3. 20) and (3.2.8) we have: 



/■ ' s 

H, 



H^ = 


= H 



^ ^‘3 , 




Using the relations 



^ Z LU^ 


— ^-1 — ■ 
and j/3 = /) cK 
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(3.3.21} can be written in UJi 


and as 


UJ^ = H ( -c^o eC, -f o<2. 0(3 ^ 

= aH 4 oCj.oCj ) 

Also satisfy the constraint 

^o<^3. a<i^ r= f (3.3.23) 


As equations (3.3.22) and (3.3.23) are linear in 



0^2 


and 

<^3 r they 

can be solved 

of 



and 

H 

, i.e. 


oCo 

C>d| 



Ujf 

-f A ) 




_ 1 
Z 

(^3 

UJy 

H 

- + 0 

(3.3.24) 

<1 

OCj 

1 

i. 

ilLl 

H 


- 1 ) 



^5 

__ / 
”” Zhi 



-f- Xy U>!^ ) 



TOG3NM1 PAGE IS 
OF POOR QUAIW 
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The ratio 


and 


can be easily calculated ; 


H — -ly 

A X, (Jj 


H 4 


— H — 

4- A Xj U/S 


_ f 


H - X^ujj, 


Substituting (3.3.25) and (3.3.26) into equation 
(3.3.18), we have 






‘ -2-T4>-/^ ^ fe-^v) ^<4 

hi 4 ^'y ^ 4 3^ 

' ~ i ~ 4- N (^g - -Ty) 


I H- Xj 




H - icy 




_L r^; ^T4N M;y ^ 

-2-L H-+:^t<y H-4j:^icy ^ ^ 

^ L H - Xy uy H - div vo^ J 3 


H dy 


Now we have four decoupled, time~ variant homogenous 
linear equations. Their solutions are immediately ' 
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available in terms of quadratures . 

Also because (3.3.27) have periodic coefficients 

(w , 0 ) , 0 ) are periodic), by Floguet theory [36], the 
X y z 

solution can be expressed in the form Q{t) exp(st), where 
Q(t) is a periodic function and s is a constant. With 
this in mind, solutions for (3.3.27) are 


oCo(^) = 
cC^ (Jt) = 
(;t) = 


Jhi-p (JiRi) 
(-Ap) (-ARt) 


where 


E, = 


H -f Jy \Jy(±) 

H Xy OJy [ tc) 




H - Jy 


H- Xy i^y 


J, 

"z 


p. =if^- 

^ L H 


±L 

1/ J L 


*£J 


OJ 


H 
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p^ = 4-i 

< H 


Jj 

TT(TZo)-t 


R= — 

f ^ 

+ - 


> 

rr(Toj) 


2, 


i H 

Hy J 

~^oJ , 



and 


TT(*) =: <it: 

Jo 1 - 


b <j> 


n(t) is an elliptic integral of the third kind. 

is the period of u , m , and w • 4> and t are re- 

X y z 

lated by equation (3.3.6). Also ’('(t) is given by 







Note there are two frequencies involved in the 
solution of Euler symmetric parameters : the first one 
is the same as tlie angular velocity and 

with period of T<n , and the second is related to exp{iA,R>() 
with period of . The latter one can be explained 

as due to the motion of the axis of the instantaneous 
angular velocity vector W,, . In Poinsot 
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vector makes a locus on 


construction, the tip of the 
the invariant plane called the herpolhode ^nd also a 
locus on the momentum ellipsoid called polhode. The time 
required for the vector to complete a closed 

polhode locus is . We call *R’ the polhode frequency. 

Prom ^ r the general solution for 

is : 


✓ \ 


V 

EjCc^CR+Ri) , o ^ -6|^(Pi+Rt}, ° 

p.{^> 


o ^ EiCcURi+R*) , ® . Ei-At^CPi+Rt) 



E, Al»vCp,+Rt}, o / &, o 



0 , -e^^tpa+Rt) , 0 , Ei'^CPz+Rt) 


PJ-ib) 


( 3 . 3 . 29 ) 


X- 

(C) The Transformation Matrix Cub 

'X 

The transformation matrix C^l, in terms of 
is 
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where 3i is given by (3.3,291, BY direct substitution, vre 


have 


a 


^ f / 

Clj + Cj, 


(3.3.30) 


where 




^ -2 0=0 2f, Bi"' ^J<^) 


•5f/ 

C, 


.ie.E^f C«r,Cp;+p^;(^,(o) 

K(o)^^A^) h^o)) 

•+^(h+P:i)[|^^(c;/3^Co) 


0i(o)-l^, f ‘>>1) 


-J ^15. if f'ii‘jPih) 
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o* 


/ 


— 2 p, E,^ ( )) 

-2 ^(cj ^Jc) ^ 

^ ^1 ^2 f-M^Cpi-f-Pj) { 

+ c<T5 (5+p^) ^ 0 

f /3ifc)i3j^<?)4 A^o)^Jo;| 

^2 [-2 >4nx.2{^ f^i(*>)P^to) 

4 C«o 2 p^ I' ^3^ f 0} -y3,^(ojj) 

O 

Ej * 

•(■ C^ 2p, 


Ei|^2^2pj, p^^o)(i^(cy 

H- ^P2 f 

^ •+ E,""{ 2^2 lp, fiU^)jiJoJ 

-2 E, E 2 LPt^F0(/z(eJj!^cJ 
) + ^Cp|4Pt) 

i^cco t 

e,’’ [ 

“”2Cc«i2p, E,^ 

02^ C/33V)7^>)] 

— 2 ^iCo)pj(c) 


^ ZE,E^ Cp,-P^) ( filo)/^^o)-^^,lo}^,co)] 

+ Ai^ ^?^-p^,)[|}^l„)^,(o)-pJo) I 

0 , E,^(i?^o)4p^^D)] 

^ ^ ^ 6i E^, [ C<J^ Cp,-pj A^W. 
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t * t 

Note that matrices C 2 and are periodic functions 

with period of only. 

Summary of the Section 

1. Without external torque, both Euler's equations and 
Euler symmetric parameters can be solved analytically. 
The solutions are given by equations (3,3.7) to 
(3.3.10) and (3.3.29). 

2. The angular velocity is a periodic function, 

that is Tio)^ , whereas P i±) 

contains two different frequencies; they are the w — 
frequency and the polhode-frequency. 

3. The transformation matrix is given by equation 

(3.3.30), in which the ^ -frequency and the 
polhode-frequency are factored. 




OF 






IS 
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3.4 Disturbing Torques on a Satellite 

Of all the possible disturbing torques which act 
on a satellite, the gravity gradient torque (G'.G.T.) 
and the geomagnetic torque (G.M.T.) are, by far, the 
most important. Fig. 3,4.1 illustrates the order of 
magnitude of various passive external torques on a typical 
satellite [17] , in which torques are plotted in terms of 
altitude. Note that, except for very low orbits, the 
G.G.T. and the G.M.T. are at least a hundred times as 
big as the others . 



Figure 3.4.1. 

Tor<iues on a sateUite of the Eaith as a function of 
the orbii height A: gravity torque; Af* aerodynamic torque, 

magnetic torquejitfc solai radiation torque, 
miccomeieoriie impact torque. 
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(A) Gravity Gradient Torque On A Satellite 


In this section, we first derive the equation for 
the G.G.T.; second its order of magnitude is discussed; 
third, by re-arranging terms, we express the G.G.T. equation 
in a particular form by grouping the orbit-influenced 
terms and the attitude-influenced terms separately, so 
that it can be handily applied in our asymptotic analysis. 

1. Equation for Gravity Gradient Torque 

with ^rf ° 3.S defined in Pig. 3.4.1, 

let R == |Rl f and further define f = R , and r == Irl . 



Hi 


X; 



C- center of mass 
R- position vector 


Fig. 3.4.1 
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The gravity attraction force acting on the mass dm by 
the earth is : 

dF = -u dm r/r^ (3.4.1) 

where u is called the gravitation constant (u = 

1 fi 9 

1.407639 10 ft /sec ) . The corresponding torque generated 
by dP with respect to center of mass c will be: 

ctL = f X 

—AX c( )n / — ~ 

« — p— (f xr) 

= [ f x( R + 

= (f^R) (3.4.2) 

We note that: 





r* r ^ ) 

^ R J 




1 


+ — 


(3.4.3) 


We have omitted all the terms smaller than 1/I& . 
The total torque acting upon the satellite using (3.4.2) 
and (3,4.3) will be: 
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By definition, because c is the center of mass, the 
first term in the above equation has to be zero, therefore: 

J ( 3 . 4 . 4 ) 

If L is expressed in body principal axes, let: 



Then, 



3 /U 

~w 


f r K,ejr, ■ 

W.A + R^eP, + RaR,f/ 
Z-R4R1S.P4 Jl. ?>££>. _ 
''iRifr + + R.iy;£ 


Using the fact that, in body fixed principal axes, all 
the cross product moments of inertia are zero, 

=0 it 

We have, therefore, 

C ) 

T b ^ ^ ^ 

l-% ^ Ri (-3.4.5) 

_R, - X^) 
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In vector notation : 


where 



(3.4.6) 



is the matrix of 
moment of inertia 


and 

^ (fx- 

^/ = I ( f/- f,") 

f/) 

(3.4.6) is the equation for the gravity gradient torque; 
it has a simple format. 


2) Order of magnitude considerations 

Since G.G.T. is the major disturbing torque for satellite 
attitude motion, it is important to know its order of magni- 
tude. As a matter of fact, this information - order of mag- 
nitude - plays an essential role in our analysis of the dy- 
namics by asymptotic methods . 

From orbital dynamics [37] , the magnitude of the posi- 
tion vector R, can be expressed in terms of eccentricity e 
and true anomaly f, that is: 

R = - y — , (‘3.4.7) 

1 -t fe. -f 
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where 'a' is the orbit semi-major axis. The orbit period 


•p' is; 



Z,7T 


u>. 


CrhtX 



(3.4.8) 


Combining (3. 4. -6), (3.4.7), and (3.4.8) we have: 



3 CO, 


CfH.t 


(/•# 


(I- ^"3" 



O ( } 


(3.4.9) 


We can say that the G.G.T. has the order 

of the square of the orbital frequency, if the eccentri- 

city e is far from one (parabolic if e=l) , and if the moment 

* 

of inertia matrix is not approximately an identity matrix 
{special mass distribution) . 


3) Re-grouping Terms for G.G.T. 

We find that - the orbital position vector expres- 

sed in the body fixed coordinates - contains both the orbital 
and the attitude modes. Inspecting the G.G.T. equation 
(3.4.6) , it seems possible to group the orbital and the atti- 
tude modes separately. 

Since 



R. 


) 

^ -/d*rv4 / 


where 'i' denotes perigee coordinated, we write 
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— k * ^ 


7 (I-*- €<^-f ) 


I CtA-^ 

t.J ° I" 


* ^ f^-f 

■^tn I O 

(3.4.10) 


Let us define a linear operator OP(B) on a matrix B by: 


OP C6) = 



o 

6z, 

B, j 



o 

Ba, 

^Bi3 


3 ^ 


(3.4.11) 


It is straightforward to show that: 



jc^ i] 


f 


o 

- 

X 

B 


) 

. 

\ 



f f \ 


= Of>(&) 


/Wl-J- CoQ,-f / 


(3.4.12) 


Substituting (3.4.12) into (3.4.10) we have: 
a V. * ^ X 




. -3 ‘^via(l+ec«+/c-^''f 

^ 


'7 


(3.4.13) 


From the above equation it is clear that the first 
group contains terms influenced by the attitude motion with 
higher frequency, and the second group is influenced by the 
orbital motion with relatively lower frequency. 
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(B) Geomagnetic Torque 

A satellite in orbit around the earth interacts 
with the geomagnetic field and the torque produced by 
this interaction can be defined as a vector product 




(3.4.X4) 


where B is the geomagnetic field and Vj^ is the 

magnetic moment of the spacecraft. The latter could 
arise from any current-carrying devices in the satellite 
payload as well as the eddy currents in the metal structure, 
which cause undesirable disturbing torques. On the other 
hand, the vehicle magnetic moment could also be generated 
purposely by passing an electric current through an 
onboard coil to create a torque for attitute control. 


If the geomagnetic field B is modeled as a dipole, 
it has the form [38] . 

^ Ms 


6 


/? 


I R^eg-3 (3.4.X5) 


where 


is a unit vector in the direction of the 


geomagnetic dipole axis, which inclines about 11.5 
degrees from the geophysical polar axis. Vector R 
represents the satellite position vector, is the 
magnetic constant of the earth ( Mg = 8.1 ><10 gauss-cm^ ) . 
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Combining equations (3.4.14) and (3.4.15) and expressing 
(_i^ in body fixed coordinates , we have 



(3.4.16) 

Although neither the geomagnetic field nor 
the body magnetic moment can be determined precisely 
in general , modeling both of them as 
dipoles will be sufficiently accurate for our purpose. 

Summary of the Section 

1. Gravity gradient torque (G.G.T.) and geomagnetic 
torque (G.M.T.) are by far the most influential disturbing 
torques on satellite attitude motion. 


2. The basic equation for G.G.T. is: 



(3,4.6) 


For G.li.T. is 






(3.4.16) 
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3. G.G.T. and G.M.T. have the same order 

3. 

as ^orbit r if the eccentricity is not too high 

and if the satellite mass distribution is not too nearly 
spherical . 

4. By re-grouping, we separated terms of attitude 

. r ^ T h 

frequency and terms of orbital frequency in L ^ and Lj^ 
the results are (3.4.13) and (3.4.16). 
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3 . 5 Asymptotic Approach for Attitude Motion with Small 
Disturbances 
(A) Euler's equations 

Euler ' s equations in vector form are ‘ 

W + C wx J UJ = l^Tj i- H (3.5.1) 

assuming the initial condition is 
W (:tc) == 

where Tx + ‘ ' represent the dis- 

turbing torques. The order of magnitude of these dis- 
turbing torques are discussed in the previous section. The 
small parameter E is defined as the ratio of orbital 
and attitude frequencies. 

Let be the torque-free Kirchhoff's 

solution which satisfies the particular initial condition J 
that is 

^ ^ 

ITjryx (AJ^ -h [ X j -^vn. ~ 0 (3.5.2) 

By Encke's approach [37]/ let 

lU(;(r) = -h L (3.5.3) 
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Substituting (3.5,3) into (3.5.1) and subtracting 

(3.5.2) , we have the equation for 

it ^ — ^ ^ 

+ I ^Nl ^ J Jm. 

+ £. [ $10X, ] Ijvv t Tj + t^Tz + “ • • 

(3.5.4) 

W = 0 

Note that Encke ' s perturbational approach is not an 
approximate method. Because by combining (3.5.4) and 

(3.5.2) the original equation can be reconstructed. 
Nevertheless, performing the computation in this perturba- 
tional form, one has the advantage of reducing the numeri- 
cal round-off errors . 

For simplifying the notation, a periodic matrix 

-if 

operator A W with period of Tuj can be defined as : 

f\ " — i>v\. L C C ^ 3 J (3.5.5) 

Eq. (3.5,4) can be re-written as 

JL * -1- 

%w -- A(^) + 1 I^W (^WJ.X) Xy^ $u) 

^ _i 

= 1 I,,. T| -f 1x^ X2-+ 

(3.5,6) 
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We see that the above equation is a weakly non-linear 

equation, because the non-linear terms are at least one 

order smaller than the linear ones. 

Further, by Floquet theory , a linear equation with 

* 

periodic coefficient, such as in (3. 5. 6), can 

be reduced to a constant coefficient equation, as 
follows. 

* ^ 

Let matrices and (■^) be defined as 

K. I 0)] 




where 


(3.5.7) 

JC ^ 

is the transition matrix for , 


It can be proven that [36] ^ 


( 1 ) 




(t) 


is a periodic matrix 


Pa"' (* + Tw) = Pa''(*) 


( 2 ) 


Pa(*) A(*) Pa"'(*) 4 4^' 




a constant matrix 


64 



Let Sco(^) be transformed into by 

— — — 

SUJ (3.5.8) 

Then 

m 

^ -if- ~ .Z. 

^ ^ Pfy + Pa 

“■ ^ P/S) + X/h ( i X,^ SoJ 

+ t i;,' T| + % +-- - J 

~ (Pa Pa -^ Fa a Pa ) 

^ ^ ^ _, * , ^ _1 — V 

^ U X j ( 1^ 

^ FaXw. T, + £^Fa X^ T^+“- • 

^ ^ ^ -I / * -1 - \ 4 / o"' n \ 

= Pfl U. - £ P^ ij (Pa Pa 

■)( ^ ^ — I — . 

-j- £ Pa Xy^' T, + I'^Pa Ta + '* - (3.5.9) 

These results can be proved by using eq. (3.5.5) and 

the property (2) . 

Moreover, in our case, the constant matrix 

is proportional to J/t«;> therefore it can be considered 

■K" 

to be of order £, . Also, for simplicity, 
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be transfornied into a diagonal (or Jordan) form 


Such as : 


M M 


■K 

A 


o 

o 


'X. O 


O 


3j 


(3 


Let 


ir=: M ' u. 


(3 


Substitute V~ into (3.5,9); finally we have 


vr= lA £ 0.(X^' e^'v 

+ a (£T| + £?Ti+-"- ) 


(3 


where 


* . « 


-1 


(3 


.5.10) 


.5.11) 


. 5 . 12 ) 


.5.13) 
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MTS Asymptotic Solution 


For solving (3.5.12), we apply the linear multiple 
time scales asymptotic technique [3] . We first expand 
the time domain ii into a multi-dimensional space ^ 


t X, , ) 

and 'To t » ...etc. are defined as ' 

OTo = A ; == • aU, 

The time derivatives in the new dimensions are 
transformed according to 


— ^ ^ I . ^ 1 J 

a To 


5 

^^0 


+ 1 



1^' 


5-c, 


(3.5.14) 


Also we assume that the dependent variable 
can be expanded into an asymptotic series of ^ ‘ 

y {:(•)=: X, '- ) + £ 1'^ ('To, T,) + ■ (3.5.15) 
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Substituting (3,5.14) and (3.5,15) into (3,5,12) 
"we have ‘ 


3 Tc I 3 TT/ a To / ' 3 


+ lii . ^ 


) 


-f. 


^ r ^ ^ -jf ^ ^-1 _ ~ -if-- , 

- ^ ^ - R (x^ Q. Q. +^77] 

,,. x _ ^<^1 _ 

[ A a u;=')6l V„ 

~ G. (Xm Q + ^ J'*'*” (3,5.16) 


By equating the coefficients of like powers of £. on 
both sides of (3.5.16), we have: 


Coefficient of 

— o 

a -To 

Coefficient of 
3 To 


.0 , 


(3.5.17) 




■?< — ^ ^ •‘f-l - T ~ 

+-A iXo -a[jJ fi. lt. 


a-c 


•+aT, 


(3,5.18) 
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-X . 

Coefficient of c. , 


^ T(j 




■+ (i,:' a' If X) a' C>. J + ^ (3..5.19) 


Solving the partial differential equation (3,5.17), 
we have 

=r ^oC'Tl ) (3.5.20) 

where is not a function of •7:’^^ / and the initial condi- 
tion is 


\Tc (o) = O 


and ) is to be determined. 

Note, (3.5.20) implies that the zeroth order solution 
for Euler’s equation with small disturbances is the 
Kirchhoff's solution. 

Substituting (3,5.20) into (3.5.18), xve have 










<! Xx)ft v; 


+ &T, 


(3.5.21) 
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This equation with the initial condition O 

can be solved as. 



■+ 





( 3 . 5 . 22 ) 


^ \ if-/- 

The term L/3 in the above equation 

can be written, as shown in Appendix A^as follows 1 




F;(T»>^. u; 


( 3 . 5 . 23 ) 


where F I > r and Fj are periodic matrices 

with period of Tu , and tg. , jL=/,i,3 are three 

* 

components of the vector . If we expand F, , 

* 

Fi and F3 into Fourier series 
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WQ have | 


a * _i X , _ 
a a~‘ tr, X ) 


a V, 



E^' j a 4 v ( 


Tw 




(3.5.25) 


where 

Also, 




and 


a 

F- 


are constant matrices. 


in case that the external torque T, is not an 

X- 


explicit function of time, then flL T| can be transformed 
into a particular form as. 


a 






(3.5.26) 


where ■ , ^ 2 . functions of slow variables 

X- 

only, and Qp ('tb) is a fast tiitie-vaiying 
function. The expressions of (3.5.26) for gravity gradient 
torque and geomagnetic torque are given in the following 
sections. 

Substituting (3.5,25) and (3.5.26) into (3,5,22), 




_ f 

^ r 


X = J 




r|, 
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(3.5.27) 


Vo. iTo 4 ^2 t'^i ) dTo 


For achieving a uniformly valid approximate solution, 
it is required that II i/} be bounded 

uniformly for all time . Applying this condition to 
(3.5.27), it requires the first bracket to the right 
of the equal sign to be zero since, otherwise, it will 

r 

linearly increase with Xj, , This leads to the equation: 


•avo 


A 


5 _ _ 

-I FXo ^,^■^1) 

(3.5.28) 


tut-ch JC, c. 


u; ro; = o 


Also eq. (3.5.27) becomes 


ir, (Toy-c,) ^ 

L in Ar 


(3,5.29) 
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Egs, (3,5.28) and (3.5,29) yield the first order 
multiple time scales asymptotic solution to the equation 
of f (3.5.12). Once l/'l't) is obtained, the 

asymptotic solution of Euler's equation (3.5,1) can be 

c 

constructed as ; 

(U (i:) = 

^ (3.5.30) 

^-1 ^ 

The matrices 1^ and M are given by (3.5.7) and 
(3.5.10) . 

The above approach gives us an alternative way of 
evaluating a satellite angular velocity u(t) instead of di- 
rect integration of Euler's equation. By this new approach we 
integrate (3.5,28) directly and evaluate (3.5.29) analyti- 
cally. Since (3.5,28) is in terms of , 

it allows us to use a large integration time step, 
thus saving computer time. Also, we note that (3.5.29) 
describes the oscillatory motion and (3.5,28) describes 
the secular motion of the satellite angular velocity. 
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(B) Solution For Euler Sytaitie-bric Parameters 
Euler syiTimetric parameters are kinematically related 
to the angular velocity by a linear differential equation 
(3.2.9) 



( 3 . 5 . 31 ) 


Fn.CTvv, (3-5'-3^) 

10 ~ uJtJ -h % 

= M (u;t ev^) 

Hence 

(5 = ^ ^ ( 3 . 5 . 32 ) 


Again by linear multiple time scales method 


To = ^ 
T, - tit 

JiCc 


we have 


d. _ ^ ^ £ 5 

cL± ~ 3Xo 


( 3 . 5 . 33 ) 
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Assume that 


/3oC^;'q-) + ^ /Si (%,%••)+ (3.5.34) 

Substituting (3.5.33) and (3.5.34) into (3.5.32), 
and arranging terms in the power of £, , we have ; 

Coefficient of : 






Coefficient of L 


I . 


^ 4- (3.5.36) 




(3.5.35) 


Coefficient of 


5~o ^ 




i 9- 


(3.5.37) 


StCc. 

Let 

LHvJ * that is 
K 


be the transition matrix for 






The expression for ±-p^ is given by 

(3.3.29), that is from the solution of /?,’ by Morton's 

t A— 
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approach. Similarly, M can be also achieved 

by using the Floquet theorem^ although the latter one 
is more numerically oriented and requires several trans- 
formations. 

The solution for the order equation (3.5.35) is 


Po t'Cc, ”^i ) — ("to ^ o) ("r, , Tj.- - • ) 


(3.5.38) 


with I.C. 




and ^oN determined. 

Substituting eq. (3.5.38) into (3.5.36), we have 


Ml 

at. 


* r- 




e^f 

(3.5.39) 


The solution to the above equation is : 




(3.5.40) 


where $UI {X,^Xj) , given by (3.5.30), is a function 
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of %J and 'T:, , and ^%>^) is a function of 'Vo 

only. By re-grouping terms it is possible to write 

V -K ^ 

jS^jj follows (appendix B) : 


(S'^1 % = R,(T,) + Pg(T„) R^(X,) 


Substituting {3.5.4() into (3.5.40), 


(3.5.41) 


/ 3 ,= lcTi,oJ 'C. 4 I- 

+ -kr Pg {<r> ^ 


(3.5.42) 


In order to have |l Pill/il ^ bounded in it is 
necessary that ^ j should not increase with time 
faster than /S^ . Therefore in eq. (3.5.42), those 

terms which increase linearly with -r^, have to be set to 
zero. That is 


4 ^ R| (%) ftN (-r,) = o 




(3.5.43) 


and (3.5.42) becomes 


/9j = (% , “) r p R (O ) Aw (t, ) <( cr 


(3.5.44) 
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Equations (3,5.43) and (3.5.44) give an asymptotic 

.solution to the equation of Euler symmetric parameters . 

* _ 

-c,) ^ (To / o) Cr,J + £ (To, X, ; + • 

(3.5.45) 

where from (3.5.43) gives the secular 

variation of the perturbational motion and ^ii"Cc) 
gives the non-biased oscillatory motions. 

Summary of the Section 

(A) The rigid body satellite attitude dynamics 
are described by; 


* _ 
T>v\ U) 


(3.5.1) 

9 

^ — 
r tuj /3 

(3.5.31) 


which can be replaced by the asymptotic approximate 
formulation 


^ Ji u; - X 




(3.5.28) 
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(3.5.30) 


__ ^ — I a ' 

UJi-k)=^ e M ( -t- € ) 


C*AU^ 


A, = O 


5 X 


(3.5.43) 




/?(*)= lp,tt) 


(3.5.44) 


(3.5.45) 


(B) The secular eqs . (3.5.28) and (3.5.43) have to 

be integrated in ^ t , €. ) , or equivalently, in 

a very large time step in t . The oscillatory equations 
(3.5.29) and (3.5.44) can be analytically calculated if 
the external torques Tj are not explicit functions of 
time. 

3 . 6 Attitude Motion With Gravity Gradient Torque 

The gravity gradient torque , given by (3.4.13), 

is *. 



c ^ ^ 

= Cxu op ( X.. 


vb )] ■ [ 


3 ^gyyX 

0-ey ■ 


L Cd-p _■ 


(3.6.1) 
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where . e, f and the satellite orbit 

eccentricity, true anomaly and averaged orbital angular 

if: 

velocity, respectively. Op is a linear operator given 
by eq. (3.4.11) . 

From (3.3.30) 


^ ^ / 

CeiiC^-di) + 


Cj. + 


■X / 
^3 


and 


(3.6.2) 


-K 


A 

^ h b C yj 


i( f 

c. 


*, / 
on (2-fet) -i- 


Aaa\ (^2.12 "t ^ + 


X / ^ 

Ci Q,, 


■X • 

' Cj Cc;}(2Fit) •+■ (.’2.R t ) 4 - C3 

(3.6.3) 

X ^ ■ 

where Cj , and Cj are three periodic matrices 

with period of 

Substituting (3.6.3) into (3.6.1) 
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+ Opf'Ci 




O-feccA-f) 

U-e^)^ 


/4^-f Ctrl, p , 


t>V' 


^ ~ s,(t) + >i^ (2.R-() + Sj W) Cca,(ZRt ) 

+ Sq. tt ) >wCiffrt) + 


cohere^ 


3 04 e Ofi-f j I ^ -f 




C/-e") 



(3.6.4) 











4 


Ci Ci 


X 


jfT"^ ^f\ 

4C^ IM 


•X X 


T ^ 


+i a.%.°P ( c/x 


if if 


+ i a c, <5p ( c,’’i„ c 


Vn 

if 

3 
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stc. 




Because -S^ ^ J are periodic matrices 

with period Toj , they can be expanded into Fourier 
series with period of Tu» . That is 

=|J 


{ 3 . 6 . 5 ) 


With that, & can be expressed in the form: 


T-t 


a Ct, ) -f 


( 3 . 6 . 6 ) 




3 ~ j 


C&l^ -j- 
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+ • ‘ 


(3.6.9) 

^ 

and analytically integrated. 

Resonance 

in’ 

In special cases, if the tumbling frequency -=^ and 

ttij 

the polhode frequency R are low order commensurable, 
that is, there exist two small integers n and m such 
that 




(3.6.10) 


then resonance will occur. For handling resonant situations, 
for example, if — — — 2-g — O , then in equation-(3 . 6 . 9) 
of , terms such as , cy 

etc. will produce constants (i.c. ~ 

Tuj 

zz J LcjriJ-aff-f- Ccry<^ ~ 4- ). These constants should 

X X ^ ^ 

be multiplied by and grouped into in (3.6.6). 

Then our theory, with the same formulation enables us to 
predict the attitude in the resonant case as well. 
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A Numerical Example 


A rigid body asymmetric satellite in an elliptic 
orbit is simulated. Here, we assume that the satellite 
attitude motion is influenced by the earth- gravity gradient 
torque only. 

The numerical values used in this example are the 
following: 


Satellite moment of 

inertia: 

= 39.4 

slug-ft^ 

ly =33.3 

slug-ft^ 

=10.3 

slug-ft^ 

Orbit parameters of 

the satellite 

eccentricity < 

3 = 0.16 

inclination i 

= 0 

orbital period 

= 10,000 sec 

Initial conditions i 

are : 

= 0.0246 

rad/sec 

It 

o 

rad/sec 

o 

II 

rad/sec 


Satellite starts at orbit perigee and its initial 
orientation parameters are 
= 0.7071 
/3, = 0 

=0.7071 
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In this case, the small parameter £ is about 


^arbi-t 
^ att/tade 


0,03 


With the above initial conditions, the satellite attitude 
dynamics are first directly integrated using a fourth 
order Ruhge-Kutta method with a small integration time 
step size of 10 sec for a total interval of 8000 sec. This 
result is considered to be extremely accurate and referred 
to as the reference case from here on. The other simulations 
are then compared to this reference case for 
checking their accuracy. A number of runs have been 
tried, both by asymptotic approach and direct integration, 
with different time step sizes. They are summarized in 
Table 3.5.1. The errors in each case - i.e., the differences 
between each simulation and the reference case, - are 
plotted in terms of time and given by Pig. 3.6.1 through 
Pig. 3.6.12. Fig. 3.6.13 is a plot of the maximum 
numerical errors as functions of the step size. From 
this plot, we see that with direct integration the 
step size at should be no greater than 25 sec. On the 
other hand, for asymptotic simulation, the step size can 
be as large as 500 sec although the first order asymptotic 
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approximation errors approach zero as £ -> 0 but not 
as 0. Fig. 3.6.14 is a plot of the required 

computer time in terms of the step size at • We note 
that for extrapolating a single step, the asymptotic 
approach requires about double the computer time using 
direct simulation. However, since the former 
allows use of a large time step, overall, this new approach 
will have a significant numerical advantage over 
direct simulation. In our particular case, the saving 
is of order 10. Although, in this comparison, we did not 
include the computer time required for initializing an 
asymptotic approach by calculating Kirchhoff's 
solution and some Fourier series expansions etc.^ we 
argue that this fixed amount of computer time 

required for initializing (about 40 sec for 

the above example) will become only a small fraction of 
the total, if the prediction task is long. For example, with 
the. above data, if we predict the satellite attitude motion 
for an interval of three days, the direct integration 
with a step size at = 20 sec requires 1700 sec of computer 
time, while the new approach with Ay = 500 sec needs 

about 170 sec plus the initialization of 40 sec. 
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3.7 Satellite Attitude Motion With Geomagnetic Torque 
The geomagnetic torque acting on a satellite can 
be written as, (3.4.16) 





^8 



(3.7.1) 


b • » 

where is the vehicle magnetic moment and 

Gjj is the geomagnetic dipole axis. The latter , 
for simplicity, is assumed to be co-axial v;ith the 
earth north pole. 

Substituting the expression of CL;Lb from 
” t) ^ — b 

(3.3.30) into Lj^ , we find Q. can be written 

in the form of 


S- = S-p (IT.) Sj. (.-r, ) 


(3.7.2) 


where 

X ^ • b I ^ 

6Lp Ct,) =CK Vm ^ 

K 

(3.7.3) 
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r>*. 



i f ^ i 

By expanding Cj , C^' , into Fourier series, 

X -k 

■■ we see that (i can be analytically integrated in terms 
of -To . Eq. (3.7.2) corresponds to (3.5.26) in section 3.5 
and the asymptotic formulations can be easily applied. 

Numerical Example 

A rigid body satellite perturbed by the geomagnetic 
torque is simulated with the same satellite which flies in 
in the same orbit as given in section 3.6. In addition, 
we suppose that the vehicle carries a magnetic dipole , 
which is aligned with the body x-axis^ 

'<Xmp Se,c 

At mean time, the value of the geomagnetic 
field is assumed to be | 

See - 

Using the above numbers, we simulated the satellite 
dynamics by direct integration and by the asymptotic 
approach. Table 3.7.1 lists all the runs tried. The 
errors of each case are plotted in Fig. 3.7,1 through 
Fig, 3.7.12. Similar conclusions as given in the gravity 
gradient case can also be reached for the case of geomag- 
netic torque. 
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CHAPTER 4 


PREDICTION OF ATTITUDE MOTION FOR A CLASS OF 
DUAL SPIN SATELLITE 


4.1 Introduction 

In chapter 3, a technique for speeding up the pre- 
diction of a rigid body satellite attitude motion was 
developed. However, the limitation that requires the 
satellite to be a rigid body seems severe, because many 
satellites in operation today have one or several high 
speed fly-wheels mounted onboard for the control or 
stabilization of their attitude motion. The combination 
of the vehicle and its flywheels sometimes is referred to 
as a dual spin satellite. Therefore, it seems desirable 
to expand our prediction method for handling the dual 
spin case as well. In doing so, it turns out that— it 
is not difficult to modify our formulations to include 
the dual spin satellite, if the following conditions 
hold : 

1. The angular velocity wJ is a periodic function 
when there is no external torque, 

2. A torque-free analytic solution of the system 
is possible. 

3. External torques are small. ' 

However, the dynamic characteristics of a dual spin 
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satellite in general are not clear yet as far as conditions 
one and two are concerned. Although we believe condition 
two might be relaxed by some further research, they are 
beyond the scope of this effort. 

As a demonstration example -for handling a dual spin 
case, we consider a special class of dual spin satellites: 
that is, a vehicle having a single fly-wheel which is 
mounted along one of the vehicle body principal axes. 

The satellite is allowed to have an arbitrary initial 
condition and to move in an arbitrary elliptic orbit. 

In what follows, the rotational dynamics of a dual 
spin body are first discussed. An excellent reference 
on this subject is by Leimanis [22] . Later, the torque 
free solution - which serves as a nominal trajectory - 
for a class of dual spin satellites is presented. This 
torque-free solution was first given by Leipholz [23] 
in his study of the attitude motion of an airplane with 
a single rotary engine. Then, in section 4. 4, an asymptotic 
formulation for a dual spin satellite is discussed and 
two sets of numerical simulations are presented. 
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4.2 Rotational Dynamics of a Dual Spin Satellite 


For a dual spin satellite, we assume that the 
relative motion between fly-wheels and the vehicle do 
not alter the overall mass distribution of the combination. 
Thus, for convenience, the total, angular momentum of 
the system about its center of mass can be resolved into 
two components. They are: H , the angular momentum 

due to the rotational motion of the v/hole system regarded 
as a rigid body, and the angular momentum of the 

fly-wheels with respect to the satellite. 

The rotational motion of a dual spin satellite 
is, therefore, described by: 


dt 


(H + H„) = 





10 -h = 


(4.2.1) 


* 

where is the moment of inertia of the combina- 

tion (vehicle and wheels) , and is the external 

disturbing torque. 

Applying Coriolis law, we can transfer the above equation 
into body fixed coordinates 'b' , 


it 

a 




d ± 


■ — Aj _ 

■+■ aj X Xi 




U) -h 




rr M 


(4.2.2) 
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Further, by assuming that the wheels have constant 


angular velocities with respect to the vehicle. 




O 


Hence 




(-IbA ) = M 


( 4 . 2 . 3 ) 


This equation is equivalent to Euler's equation 
in a rigid body case. 

For satellite orientations, because Euler syimnetric 
parameters and the angular velocity ai are related 
kinematically, the equation remains the same for a dual 
spin satellite, that is. 



4 . 3 The Torque- free Solution 

The motion of a dual spin satellite with a single 
fly-wheel mounted along one of the body principal axes 
can be' analytically solved if there is no external torque 
acting on the vehicle. 
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(A) Solution to Euler’s equations 

Without loss of generality , suppose the fly- 
wheel is aligned with the principal x-axis . By 
(4.2.3), the Euler's equations are 


c- 


where 'h' is the angular momentum of the wheel with 
respect to the vehicle. 

If eq. (4.3,2) is divided by (4.3.1), we have 

Integrating this equation and solving for uiy in 
terms of ^ the result is \ 



( Ij - ij( ) ±y 10^ —2. Xy lO^ 


+ C, 




(4.3.4) 


where .Cy is a constant which is yet to be de- 
termined , 
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can 


Similarly, using eq. (4.3.3) and (4.3.1), 
be also expressed in teinns of ^ "the result is 




( Jx - 2, Tx ^ Cj 


(4.3.5) 


where is another constant. 

In case that there is no external torque, we know 
that the rotational kinetic energy T and the total 
angular momentum of the system must remain constant. 

i .e . , 


j 


W/jT ~t 1} IaJ 


7 




-tI to 






(4.3.6) 

and 

-t Xy 

(4.3.7) 

Equations (4.3.6) and (4,3.7) can be derived by integra- 
ting [ 

( 4‘6= O ■+ ■+ 

respectively. 
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Having determined the angular momentum and kinetic 
energy T, one can calculate the constants c^ and c^ 
by substituting (4.3.4) and (4.3.5) into (4.3.6) and 
( 4 , 3 . 7 ) j they are f 


Cy = - -A, ~ T- 7$. ) 

^ , (4.3.8) 

=: Hy - (:iT- ^ 


or and 


can be rewritten as ' 


VO. = ± 






(4.3.9) 


Ui. = + 




+ I ,..^y p 


where 








•X, 


‘^y 


X. ^ 


(4.3,10) 


To find 
and 


one can use (4.3.9), the solutions 
with (4.3.1) to eliminate tOu and 
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/ such that an equation for can be obtained. 


AJik. =±jr^ 


Or 


7t: = 


t'Jx' cl OJx 


Suppose that ~ iy has roots of 

and in a descending order, then 


(4.3.11) 


(4.3.12) 


> ^3 


± 



ci COx 


S 





f4.-l,13) 

This integral can be easily transformed into an 
elliptical integral; i.e. eq. (4.3.13) is equivalent to 


t ^ yr\. 


4 


£>14* 




(4.3.14) 


where 


J (Jy- X;,) [X^ -X^; ( UJ, ) UJ^ ) 


(4.3.15) 



OJ^ 

L\J^ — (aJ^ 


to 2 _ ^ 

Wj ~ W- 


(4.3.16) 
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and if 


W4 < <• 0^3 f then 


VO ^ ^ ^'1 )’tUj,(uJtf. 


(4.3.17) 


and if 


^ ^\ , then 


lU^ = 


UJ2_iVJ^~l0t ) - tU^ (^4- UJj ) My^<4> 
lO^ — Cu)^-lO^) 


(4.3.18) 


(B) Euler Syinmetric Parameters 

Since the equations of Euler symmetric parameters 
remain unchanged whether a space vehicle has fly-wheels 
or not, 

f ^ ^ 'f (4.3.19) 


Once the angular velocity u) is obtained, this 
equation for the Euler symmetric parameters can be solved 
by a similar procedure as discussed in section 3.3. The 
result can be summarized in follows. 




't,C^(P,+Rt> 0 -E,A^(P, 4 R*) 0 

f^v(o) 



0 EjC^(p^+iRt) 0 Ea.yl^(Pi+R+) 


^( 7 t) 


E, 0 Coq CPi+R*^ 0 

(V®) 



, c. -E^A;..(P:. 4 Rt) 0 0 «(.P,+Rt) 

Ho') 
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(4.3.20) 





where 


E. 


H-r + Xy ^yd) 
H-J- -^Xy LOy(O) 


> 2 . 


H-p ~ ^y(^) 

H-j.~ Xy l'Jy(O) 


_j_ 

z. 


pj = O.^ 7 r(/t) - R * 

= 0 ,s Rt 

_ TT (Tuj> 

^ T(aJ 


and 

■“Tf (^-jt ) =• f . (-^T — N-p- (^x ^ ^ 

yo H-J- •+ Xy UJy 

"llr >r: r ^ ^ X ~f~~ Hy -f ^ ^ 

^ }~i-y iVy 
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TT (*) and (±) are two elliptic inte- 
grals, T and Hy are the kinetic energy and angular 
momentum of the system, and they are given hy equations 
(4.3.6) and (4.3,7). ’ is the period of 

the angular velocity w , and h is the angular momen- 
tum of the fly-wheel with respect to the vehicle, 

4.4, Asymptotic Solution and Numerical Results 

Once the torque-free nominal solution for a dual 
spin satellite is obtained, the basic procedures of the 
asymptotic approach, described in section 3.5 for a 
rigid body satellite, can also be applied to a dual 
spin case. In order to include the gyro effect due to 
the fly-wheel, a few equations in section 3.5 have to 
be changed . 

Equation (3.5.1) has to be replaced by 

^ ^ [uJXl { Hfcc] == t^T| -*■ 

(4.4.1) 

where is the angular momentum of the flywheels 

with .respect to the vehicle. 

Equation (3.5.2) is replaced by 

Xkvv •+ - 0 (4.4.2) 
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and equation (3,5,5) has to be changed . 

^ ~ lyy^ (4.4.3) 

Of course, and the transition matrix 

arise from (4,3.9), (4.3.17) and (4.3.20). 

They are the reference trajectory for a dual spin case. 

Numerical Simulations 

For demonstrating the correctness and accuracy of 
predicting the attitude motion in the case of a dual 
spin satellite, we select a numerical example with the 
following data? 

Satellite moment of inertia 



. = 30 

slug-f t^ 


= 25 

slug-ft^ 


== 16 

slug-ft^ 


A fly-wheel is mounted along the body fixed x-axis 
of the vehicle with an angular momentum h with respect 
to the vehicle , 

h =0,2 slug-fi^/sec 
We assvune that the satellite is in an elliptic ic 
orbit with 

eccentricity e = 0,16 

inclination i = 0 

orbital period = 10,000 sec 
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The initial conditions are: 


Angular velocity uJ io) 


— 0 « 0 3 

Ujy = 0.01 

= 0.001 

Euler symmetric 


rad/ sec 
rad/sec 
rad/sec 
parameters 


/3^ = 0.7071 
|S, = 0.1031 

= 0,1065 
= 0.6913 

For these numbers, the small parameter £. of the 
problem, defined as the ratio of orbital and attitude 
frequencies , is about 




orbital frequency 
attitude frequency 


• 2 -tt /io,ooo 




This dual spin satellite is first assumed to be 
disturbed by the gravity gradient torque only. The 
dynamics are simulated both by direct integration and the 
asymptotic approach. The results are summarized in 
Table 4.4,1. Also the simulation errors are presented 
in Fig, 4,4.1 through 4.4,10; they are given in the same 
way as in section 3.5 for a rigid body satellite. 
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Next, the dual spin satellite is assumed to be 
disturbed by the geomagnetic torque only, A magnetic 
dipole is placed onboard, with a strength of ^ = 
(3,0,0) f t - amp '■ sec, which interacts with the earth's 
magnetic field of strength 






Similarly, the attitude dynamics are simulated by 
direct nvunerical integration and the asymptotic approach 
The results are summarized in Table 4,4,2 and Pig, 4.4,11 
through Pig, 4,4.20. 

Conclusion 

These two sets of data, one for gravity gradient 
torque the other for geomagnetic torque, show that 
our asymptotic approach is equally useful for a dual 
spin satellite as for a rigid body case, if the conditions 
listed in section 4,1 can be satisfied. The niimerical 
advantage of saving computer time and the approximation 
error introduced by our asymptotic approach are of 
similar character as discussed for a rigid body case. 

The details are not repeated again. 
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CHAPTER 5 


DESIGN OF A MAGNETIC ATTITUDE CONTROL SYSTEM 
USING MTS METHOD 


5 . 1 I ntr oduc tion 

The interaction between the satellite body magnetic 
moment with the geomagnetic field produces a torque on 
the satellite. This torque, however, can be harnessed 
as a control force for the vehicle attitude motion. 

By installing one or several current-carrying coils 
onboard, it is possible to generate an adjustable 
magnetic moment inside the vehicle and thus a control 
torque for the satellite. This magnetic attitude control 
device, using only the vehicle-enviroment interaction, 
needs no fuel and has no moving parts, it may conceivably 
increase the reliabiltiy of a satellite. In recent years, 
it has received considerable attention. 

To design such a system, nevertheless, is difficult, 
because the control torque is very small. Since the elec- 
tric currents available to feed through the onboard coils 
are limited, the magnetic torque generated in : 
this way is not large enough to correct satellite attitude 
motion in a short period of time. In fact, it is realized 
that one has to depend on the long term accumulating 
control effort of the geomagnetic interaction to bring 
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the vehicle into a desired orientation. For this reason, 
the system can not be studied readily by classic control 
design techniques . 

However, this problem can be more efficiently 
analyzed in terms of the slow variational equation 
from the I-ITS approach. By casting the dynamics of the 
above system into an MTS formulation, the fast motion 
(attitude nutational oscillation) and the slow motion 
(amplitude variation of the nutation) can be separated. 
Even though the effect of the magnetic torque on the 
dynamics is very difficult to observe and comprehend in 
real time t , still, using the slow secular equation 
in terms of a slow clock, the control effect on tne 
nutation amplitude change immediately becomes transparent. 
In this chapter, we will analyze a magnetic attitude 
control system for a dual spin, earth-pointing satellite. 
For more imformation, the reader may refer to the works 
of Renard(39], Wheeler [40J and Alfriend [41]. 

5,2 Problem Formulation 
The Problem 

A dual spin satellite moves in a circular orbit; 
its antenna is required to point toward the center of 
the earth. A momentum wheel is assumed to be mounted 
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along the satellite pitch axis for control of the 
pitch motion. 



For the above satellite a roll-yaw oscillation, 
called nutation, is possible. This is because its 
angular momentum vector may not be perfectly aligned 
with its angular velocity vector due to external dis- 
turbance or initial misalignent etc. . 

A magnetic control device, using the geomag- 
netic interaction, is to be designed to damp out 

the nutational oscillation as well as to keep vehicle's 
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angular momentum perpendicular to the orbit plane 


Equations of Motion 

The Euler's equations for a dual spin satellite 
are given by (4.2.3), they are * 


+ ( X 3 - ly ) 4 

Xy Wy4 jp ^ U}^ 



M 


3 (5.2.1) 



where 'h' is the angular momentum of the fly-wheel 
and M is the external torque on the vehicle. 

Assuming , 4> and B are the Euler angles which 
correspond to yaw, ' roll and pitch for concatenated 
rotations, indicated by 

/trajectory 


I 


frame /dotation about X rotation about y, 

(x,y, z) {-(x. ,y ,z. ) 

by amount of V' by amount of 4 

/ body \ 

rotation about 2 ^ \ frame] 

(x-,y ,z^) 

by amount of O' b )> b 


Then, a kinematic relation between the angular velo- 
city and the derivatives of the Euler angles can be 
written as [42] : 
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(5.2.2) 



* /\ • A 


» /N A 

0 + Wc t Ej 


A A /V 

where ^ are unit vectors, is in the 

. . ^ 
direction of vehicle position vector and 0-^ is normal 

Z' 

to the orbital plane. The term ( H;t) describes 

the rotation motion of the satellite in a circular orbit. 

^ Aw 

Transforming ^ 6^ and into the body 

fixed coordinates and using small angles assximption, 
equation (5,2,2) can be written as; 


U 


b 


4^ 

zUXx + e + 


4 y' 

$■ + fAJo 


(5.2.3) 


(5.2.4) 


Substitute (5.2.4) into Euler's equation (5.2.1), 
we have ; 
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+ (ij- lu) (f 4- 


M; 

Ii G*- ~ 

^ <5 (5.2.5 

Earth's Magnetic Field 

The geomagnetic field B can be approximately 
represented by a magnetic dipole as discussed in section 

3.4. 




(5.2.6) 


where is a constant of the geomagnetic field 

( /Ug = 8.1 X 10^^ gauss-cm'^ ) » is the unit 

vector along the dipole axis and R is the vehicle 
position vector. Using Fig. 5,1, ^ can be expressed 

in the axes of c.;< z Cy ^ (called trajectory 

axes ) \ 

" >4*44^ (XJq ^ A IOq ^ A ^ y 

V ^ 

Ccro X (5.2.7) 

where'!’ is the inclination of the orbit. Substitute 
(5.2.7) into (5.2,6), we have; 


— 'tffti vWq r I 

j3 SI ■ ^ ^ ^ 




+ c«; 


(5.2.8) 


108 



Further, because the satellite is in a circular orbit 





•and if the vehicle is well controlled, the body axes of 
the satellite should be closely aligned with the trajec- 
tory axes^ therefore, 


5 ^ 


6 ^ 




or 


V 

—“2 
t<«L A 








B>£, ^ ^0^ 




{ 5 . 2 . 9 ') 




Bo =■ 




A 



AC 


Cffi. A 
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The Control Torque 

Assuming to be the magnetic moment of the 
vehicle i the torque experienced by the satellite due to 
the geomagnetic field is 

M — (5,2,10) 

For the control of the vehicle's attitude, let us pre- 
specify the control torque M as 



4 

/• > 


-At 


ly 4 > 

- ■^1 

-Ar 


O 


o 


represents a desirable control law, 
which may not be possible to implement. 

is two-dimensional, because the pitch motion is controlled 
by the momentum wheel, which is excluded in this analysis. 
The first term in the above equation^ ^ ^ O') 

reflects elimination of excessive angular momentum due 
to the perturbed vehicle body rotational motion. and 
thus damps the nutational motion. The second term 
is to control the vehicle orientation through the gyro 
effect of the momentum wheel. 

Let us take the cross product of equation (5,2.10) 
by B , 
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BxM « Bx(vi^xB) 

- B (Vm - 8 ) 


(5.2.12) 


The above equation may be satisfied if we pick 



BxM 



(5.2.13) 


In this case Vj^ will be perpendicular to B and 
also to M . It can be proved that this is 

the smallest magnetic moment required for generating a 
known torque M, Substitute (5.2.11) into (5.2.13) to 
give: 






(5.2.14) 


®x ( — hi. ^ AfAP) 

We see that is three dimensional. Thus in 

order to implement the above control law, it requires 
three electric current-carrying coils and their supple- 
mental equipment to be mounted orthogonally to each other 
in the vehicle. 

However, let us suppose that there is a weight and 
space restriction. We have to limit ourselves to use 
one coil only. By putting a single coil along the pitch 


111 



axiSj the first and second components of equation 

(5.2.14) are eliminated^ the control law now is 

r 

S 

O 
O 

' 4 

(5.2.15) 

The corresponding torque generated by this control 
law v;ill be; 




Vm X B 

ITTlt . 



(5.2.16) 


By substituting (5.2,16) and (5.2.9) into the equation 
of motion (5.2.5) we have: 


i l8l^ lil" "" 


4 -A. tdo ] 4 ^0 A 

^ 131 ^ 


s= 0 


112 



and 


I I 

XyCp + 




iBl" 


t 


*K 

' uJ, ■¥ UJ^ (Xx-^i)- it] f 

[bl" J 

' -i , , o , >4jVv2. Wj,'"^ "A I 1 it 


- - A, A f 


(5.2.17b) 


These equations are linear, but can not be exactly 
solved, since the coefficients are time- varying. The 
quantities kj and should be picked such that 

the above system is stable, 

5.3 System Analysis 

Free Response Solution 

The roll-yaw dynamics (5.2.5) are first determined 
without considering external torques. We see that 
the dynamics contain a fast nutational mode and a slow 
orbital mode. 

The system (5.2.5), without the external torque H , 
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has the characteristic equation; 






■a 










The four characteristic roots are: 




( 5 . 3 . 1 ) 




orbital mode 


•*~a 1 - nutational mode 

' Since u)o (the orbital angular velocity) is much 
smaller than H;Ast£. ( ‘ ) , the angular 

velocity of the fly wheel, therefore the nutational 


-f A 

mode can be simplified into - — 


J Ix^ 


We 


note that the frequency ratio of orbital mode and 

. , . IOq / 0 

nutational mode xs a very small number^ that xs / «• 

Order of Magnitude Consideration 
Let £. -be a small parameter which is defined as 
the ratio of orbital frequency to the nutational 
frequency. 


- 






< < 


( 5 . 3 . 2 ) 
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Assxune h, ^ t^y and are of order one, then 

the following terms are ordered as : 

A- OCf ) 

0(6) 

6c = w;/ OC^) (5.3.3) 

6 ^ =: UJcT M-q CciiX/M. = O C€X') 

{8j = o ( e-) 

The control gains kj , k^ are limited by the small 
current generating the torque and is therefore , necessari- 
ly small. We assume that ^ s O 


For convenience, 6 and its power can be used to 
indicate the order of each term. For instance Uif, 
will be replaced by . By doing so, equations (5,2,17) 


are 


+ t r ^ I I ^ + Af 


IBI 


- Xy •+{ ) cjJ^] f 

■j- ^ f ~ ^ A A tUc^+ 

L 1 ia I ^ J 


•+ 1 ( ^ ] f =o 

^ 1st ■* 


(5.3.4) 
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Bp 


Ai. 1 ^ ^ ^ 

i&r 


1 ^ J[^ uJe + ^'c(^X"-^5)J 'f 

u.,Ji-t * ^'^ -J t 


Il3j‘ 


— 4 Bp 

IBI^ 


■^1 ^ 


y =o 


(5.3.5) 


From these equations , it is easy to see that all 
the control terms are at least one order less than the 
system dynamics^' therefore the control does not influence 
the fast dynamics. 

To simplify the notation/ '(5.3,4) and (5.3.5) can be 
re-written as : 

Xx ^ -4* ^ (-^ + €- P^p ) ^ 4 E. 

(5.3.4.) 

9 -I- (- A -UPf ) r -t t E^f 

■ ^ Ff iP — (5.3.5'i 

where ....etc. are defined by 

comparing (5.3.4'^ ) to (5.3.4) and (5.3.5') to (5.3.5). 
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Multiple Time Scales Approach 

/ / 

For solving equation (5,3,4) and (5,3,5), we use 
the MTS method as outlined in section 2,2, by 
expanding the dependent variables vf ^ into 

asymptotic series in €, , 

= fc + • 

f = + -* • 

Further, the dimension of the independent variable 
t (time) is also expanded into multi-dimensions as given 
by (2.2.11), (2..2.12), and (2.2.13). Expressing equations 
(5.3,4') and (5.3.5') in new time dimensions, we have: 



-t 

:l t ^ 


+ ^ 



^ V 



ax, 




-f 

£ d . 

nr. ^ 



4 1 

+ j 



V sro 


dr, 



+ 



STj, 

-t e 

3-C, 

+ t + 

<y 1 0 

4 


C'fp+ i 


4 - 

)- 

1 ( fv-i- 


O 

(5.3.8) 


(5.3.6) 

(5.3.7) 
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and 


iy(i!5|. + 2t^^+ +....) 

^ ^Tp^T, 5%-^ ^ 

i<f ii^ ^ 1 4^ * ^^^■■■■) 

' •d TTp 9 X^ 3 Xp 

4 f-Ti+ sjv) e 4^ + 1 


a Xe 


^ C «po+ ^ •+ — ■) -f £ Rp 


= o 


(5.3,9) 


By equating 

the 

terms of like 

powers 

above equations , we 

have: 


Terms of 




3:x 

~d Tfc^ 

- + ^ 

a To 

=: 0 

i^. 


_ /) 3 y^o 

=- 0 


(5.3.10) 


The solutions - are ; 



(5.3.11) 
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but 


Inhere p, q, r, and s are not functions 'Zo , 
they may be functions of 'f, , . They are yet to 


be determined. 
Terms of 






4 -L 


9T. 


_ _ 1!^ 

>fo _ p 




^-Ci 'f' 5 To 

^ •" ff Tc 


I - 

ax,- 


-L 


ati 


* -2 X 


- (5.3.12) 








3tr 


/ 


f 2fXo 


- - F<p 


'Tp • <p 

Substituting the solutions of Sp^ and from 

(5.3.11) into (5.3.12), we have ^ 
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where 





-■L 


% 

at, 





a r 
ax, 


- Hep ^ 


K 


Note yVlf/ ... etc. are all functions of '~^i^Ti^ 
only. And p, q are yet to be determined. 

In order that \ % / %\ and j^i /<f>^ j be bounded, 
the terms that produce secular effects in 
must be eliminated. Equations (5,3,13) are 

linear^ the transition matrix ^ o J for 
^4^ and in (5.3.13) is 
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§(X. o) 


NcA-Tp — Mu.T{, 




Ce^ Nw_Ty 




-X. ( A fNJH.T.>) 


- A -i^f7(-ji Mu Tcj 

'^P(aNuX^) 

■l-J^f^(-AUuTO^ 

(5.3.15) 


The solutions for ^ with zero initial 

conditions are •; 


te' 


f-n 

^ . ■> 
Mu. <7) +-^^’fpCANuc::)Vy,+'W'^ 

S'Tb 

/ 


o 

(-A NaC^ ) "Ucii + Nu<^ 


' Jicr 
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$ (Tc, o) 


r +A - (; -^v^j 


:2- Nk ' f' :2. Nu. 

+(aK|,~-VV^p) (^^;r ^;-^^^ j Aoy>^;^j^Xo) 


N 


tc. 


N, 


+ ( + TX^ A - XC^ A- ) 'T^ 


(iXy4lxx<,>> ^p(-x;Na,;r.) ^ /y 

: 2 - Nt, ^ t/ :2 

+ tw^ 4 ;w^) 4.p(-4N i^^ _ , 

^ ' Nijl 

•^fv^4-'u:.^*-x v^ + ; 

(5.3.16) 


Notice, the terms ( 'V^^Z€ ^ a ) 
and { increase 

linearly with time, in order to have 'fj , bounded, 

these two secular terms have to be set to zero. They 


are 


V^+TXy+V^^-xXpA z =(9 


V<ji + TX,^ — A •+ A tC y/ — 0 


(5.3.17) 
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which are equivalent to 


+ A == O 


(5.3.18) 


Also, using the freedom of choosing r and s , 


we can set "Wy =• O j ^ 


E>^ = O 

“Wc^i = -|^ ■" ^ - Fcp K = O 


(5.3.19) 


■ t 

Substitute f or TX^ , "V^ f ...etc, from (5.3.14); 
equations (5,3,18) are, then; 




(5.3.20) 


y 


- F.|..li<. - A - - -£aJ-\1l =o 


(5.3.21) 


where 'p* and 'q' are complex conjugate. 

A first order MTS asymptotic solution for (5,3.4) 
and (5.3.5) thus has been constructed by using (5,3.11); 
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that is 

(-NitATo) I ( Nu.A -Zt) + V- 

(5.3.22) 

f = -^F Hu.i-r,)-+ A. 

where P/ q are given by (5.3.20), and r, s are given by 
(5.3.19). All these equations are in terms of 'X) - the slow 
time scale. In order to make and decay in time, the con- 
trol torque has to be selected such that it forces p,q,r,s 
to zero through the slow equations. The general problem in 
chapter 4 (with magnetic torque) is simpler in the present 
case because the perturbed motion is kept small by the con- 
trol. Linearization therefore, seems adequate, leading to 
analytical solution. 

Selection of Feedback control gains 

We wish to select feedback gains kj and k^ 
in system (5.3.4) and (5.3.5), such that the system is 
stable. However, because the control torque is small 
and the system is time-varying , the effect of the control 
upon the dynamics is not clear. 

The above problem was then cast into an asymptotic 
formulation. By the MTS method, an approximate solution has 
been achieved as given by (5.3.22), through which 
the design problem can be interpreted in a more convenient 
way. That is because the approximate solution (.5.3.22) 
describes a constaht frequency oscillatory motion with 
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slowly time-varying amplitude p,q and biases r, s. 
p,q,r^and s are functions of the control torque and are 
given by (5,3.19) and (5.3,20). The problem can be con- 
sidered differently in that the control torque is used 
to force the amplitude p f q and biases r, s to zero. 
In this way, since the fast nutational oscillation has 
been separated out, the selection of the feedback gains 
kf'and kj_ in the slow equations (5,3.19) and (5.3.20) 
becomes much easier. 

The solution for p (<^j ) t or q(«<| )] from (5.3.20) 


xs 


f = fp 


— Cf _ 




X 


/ 


Xy 


Jly^ 


Imaginary part 


(5.3.23) 


Thus p and q will approach zero if 

- ^ ^ ^ Q 

XXx ly 

where C<j- , etc. are given by (5.3,4') and 
(5.3.5'). By substitution, (5.3,22) becomes* 

J ly 
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or equivalently f for p and q to be decaying, we 
require 

-4, < (5.3.24) 

Also r and s are governed by (5.3,19)^ which is 
a second order linear differential equation with periodic 
coefficients. The necessary and sufficient condition 
for such a system to approach zero as Xj approaches 
infinity is that the eigenvalues of ^ ^ 'Zj) 

lie in a unit disc |xl<J where is 

the transition matrix of the system and ' T* is the period 
of the periodic coefficients. Equations (5,3,19) can 
be written as : 


f 1 


f-F^ 




ax, 




a 

a K 
‘ ax, . 



■A ■ 
> 


r 

y 





/ \ 
/I 

=: 

IBI" * iBr 






y 


IBI" ' |6l" 




= (5.3.25) 
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Note that does not depend on ^2. * 

5,4 An Example 

With the problem described in section 5.2, suppose 
we have a satellite with parameters as follows: 

Moment of inertia 


= 120 

slug - ft^ 

ly = 100 

slug - ft^ 

.150 

slug - ft^ 

Angular momentum of 

the fly-wheel (along pitch 

axis) 


h =4 

slug - ff^ / sec 

and satellite is in 

an orbit of: 

eccentricity 

e = 0 

inclination 

i = 20° 

period 

= 10,000 sec 

The small parameter 

, by (5,3,3), is then 


~ £), 0) 7 

The value of £. can give a rough indication 


as to the degree of accuracy of the approximate solution. 
Since the roll-yaw motion can be approximated by 
(5.3.22) 
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"A,- 


where p, q are given by (5.3,21) and r, s are 
given by (5.3,25). For p, q to decay, we require: 



and for r, s to be decaying, the eigenvalues 7 *.^ , 
of the transition matrix of the system (5,3.25), 

^ p must be less than unity. The eigenvalues 

and are plotted in terms of in Fig. 5,4.1 

and Fig. 5, 4. 2. We see that if 0<kj< .1, then i-t 

and Ai < t , that is r,s will be damped. 

-4 -4- 

We select kj = 4 x 10 and k^ = 8 ,x 10 ♦ 

Implementation of the control law of (5.2.15) 

with the above kj and k^. requires electric 

power of about 10 watts. 

We numerically integrated the equations (5.3,4) 

and (5,3,5) with the above kj and kj_ for an initial 

condition of 0 . 3 ® ^ =: «,/ *’ Xf&c- , 

lo) = o. . The roll-yaw motion is plotted in 

Fig. 5.4.3 and Fig, 5.4.4, With the same numerical 

values, we also simulated the approximate solution as 
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given by (5.3.22). The roll-yaw motion by this new 
approach is also plotted, in Fig. 5.4.5 and Fig. 
5.4.6. Comparing Fig. 5.4.3 to Fig. 5.4.5 and Fig. 5.4.4 
to Fig. 5.4.6, we found our asymptotic solution to be 
very accurate. Further, the system is stable in the 
way we expected and the variations of p, g and r, s 
can be clearly identified. 

5.5 Another Approach; By Generalized Multiple Scales 
(GMS) Method Using Nonlinear Clocks 
It was pointed out by Dr. Ramnath that system 
(5.3.4) and (5.3.5) can be regarded as a linear system 
with slowly time-varying coefficients. That is, the 
coefficients of the equations change at a much slower 
rate than the dynamic motion of the system [3,4]- This 
kind of problem, can be easily transformed into a 
singular perturbation problem by letting 
and changing the independent variable t into 'X. . 

Then, another approximate solution by general multiple 
time scales approach using a set of non-linear time 
scales is immediately available. 

The solution for a second order singularly perturbed 
system is reviewed in section 2.3 and a similar solution 
for an n-th order singularly perturbed system is given 
in [3,4]. To illustrate, equations (5.3.4*) and (5.3.5*) 
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can be decoupled and f” , 'Z' satisfy an equation of the form: ■ 
4- {f I 

X (tt) oC s= O 

.The GmS solution for this equation is [3,4]: 

= I tj, a, ; H.ti) 

x=o 

j . J, -k (X^nel T~ ~0 

This GMS solution employs nonlinear scales in contrast 

with the linear scales of the MTS solution. The GMS app- 
roach subsumes the MTS method and could lead to more 
accurate description of the dynamics. 

The advantages of this alternative approach are 
twofold. First, it treats the whole class of problems 
of linear slowly time-varying systems and thereby it 
can, conceivably, deal with a more complicated problem. 
Second, it is fairly easy to numerically implement 
this approximate solution which needs much less computer 
time than straight direct integration. Thus it might 
be helpful in the area of simulating a system if the 
design task has to be carried out by trial and error. 

The same numerical example as discussed in the pre- 
vious section is used. This time, the solution is approxi- 
mated by the GMS method. The result is given by Fig. 5,5.1 
for roll motion and Fig, 5,5,2 for yaw motion. The accur- 
acy is found to be excellent. Also for demonstrating that 
the new approach can save computer time, several cases have 

been tried, and the results are summarized in Table 5.5, 
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CHAPTER 6 


CONCLUSIONS AND SUGGESTIONS FOR FPTURE RESEARCH 
6.1 conclusions 

A general method has been given for fast prediction 
of a satellite attitude motion under the influence of 
its environment disturbances. An approximate solution 
to the problem is developed by using the multiple time 
scale asymptotic approach, such that the digital imple- 
mentation of these approximations would give a significant 
saving of computer time as compared to direct simulation. 
Furthermore, because the new approach has been designed 
to be very general it can handle any orbit, initial 
condition or satellite mass distribution, and so it could 
potentially become a valuable tool in satellite engineering. 

The attitude motion of a rigid body asymmetric 
satellite is first considered. By the MTS asymptotic 
technique, the slow secular and the fast oscillatory 
effects of the disturbing torque on the attitude motion 
can be immediately separated and then be evaluated indivi- 
dually. These slow and fast behaviors, combined, give 
the complete motion while divided, each describes a 
different aspect of the phenomenon. 

Similarly, a class of dual spin satellites is then 
studied, A dual spin satellite represents a vehicle 
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carrying sizable fly-wheels on-board. This model may 
resemble many satellites in use today, since the fly- 
wheel has been a common device for attitude control and 
stabilization. We have studied a special case of dual 
spin satellite with a single fly-wheel mounted along 
one of the vehicle body principal axes. However, the 
problem of a general dual spin satellite with multiple 
fly-wheels mounted in different directions seems still 
to require further research. 

The new approach is then numerically simulated 
for tv/o environment disturbances. One is a satellite 
disturbed by the gravity gradient torque and the other 
is by geomagnetic torque. The results show that the 
new method has a significant advantage over the conven- 
tional direct integration. In some situations it can 
be faster by an order of magnitude while the approxima- 
tion errors are still well bounded and acceptable. 

•A way of handling resonant situations is also dis- 
cussed. Attitude resonance will occur if the satellite 
has a mass distribution such that a low-order commensur- 
ability exists between the polhode frequency and the 
angular velocity frequency. Then there will be a sub- 
stantial increase in the secular effect due to the dis- 
turbance, We found that the resonant situation can be 


Pact? tq 
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easily detected and handled in our approach. 

Finally, as the MTS formulation separates the 
slow and fast behaviors of a satellite attitude motion, 
we use this property for the design of an attitude control 
device. In the problem, the control torque from geomag- 
netic interaction is very small. Nevertheless, its in- 
fluence on the dynamics becomes clear if we look at the 
secular behavior on a slow clock. This idea has also 
been explored in [41] , However, we believe that the 
complete solution to the problem is achieved for the 
first time and the control law is new. 

6.2 Suggestions for Future Research 

This research concerns the attitude motion of a 
satellite which is operated in a passive mode. It is 
therefore essential for us to consider all the possible 
major disturbances. We knov?, besides the gravity gra- 
dient and geomagnetic torques, there are also other 
disturbances which are important in particular situations, 
as for instance, the atmospheric drag in a low orbit 
case and solar radiation pressure for a satellite with 
large surface area. The latter is found to be important 
and worthy of research because more and more satellites, 
especially the long-lived ones, have panels to collect the 
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sun’s radiation for their energy supply. This problem, 
however, is a difficult one, since the motion of the 
earth around the sun, which changes the direction of 
the sun light, represents an even slower mode on top of 
the fast attitude rotation and the slow satellite orbital 
motion. Furthermore, the earth shadow could give a 
discontinous change of the radiation disturbing torque. 

The- second area is involved with the generalization 
of the prediction method for all dual spin satellites; 
that is, satellites which contain multiple fly-wheels. 

The difficulties of the problem have been listed in 
chapter 4. We believe that by further research, these 
difficulties could be overcome. 

In chapter 5, we have applied an asymptotic technique 
for a classic control design problem. The asymptotic 
approach seems very helpful in that particular case. 

We believe that a similar situation, whenever the control 
force is small, could also occur in other fields such 
as orbit transfer by a low thrust engine. Therefore a 
systematic study of the role of asymptotic methods in 
control theory could be very interesting and useful. 
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APPENDIX A 


EXPRESSION FOR Q(I^ ^ Q~^VX) Q~^V 


If is the moment of inertia, then for arbitrary matrix 

i( — > is 

A and vector v, the expression of (Avx) I_Av can be 


m 


re-grouped as follows: 

^ ^ r 

( A V X ) Im A V 
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'ft ic ft 

where I and OP. (A), i = 1 . 2.3 are defined by comparing the 
^ above two expressions and i = 1/2,3 are three 

components of the vecror v. With these definitions. 


a ( 3 :^ a' a ir 


^ ^ is [ o\ (1^-' f oT) 


* -f 


% (iy^o:^) ^3] ^ 


In case if Q is a periodic matrix, we can expand 

it it it 1 

Q I OP. (I ■'■) into a Fourier series* that is 




(A.3) 
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APPENDIX B 


A T A "fc 

EXPANSION OF ^>g*“ .[5oj] 3j_jj 

)( 

The transition matrix #o(t,o) has been given by (3.3.29); 
if we substitute that into the expression of , 

we will have 
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where E E- P, P„ are defined in Section 3.3. They 

i / 2 / X / 2 

are all periodic functions with period of T^. 
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Substitute 6 CO from (3.5.30), and expand all the 
^ periodic matrices by Fourier series. It is tedious but 
straightforward to have f-^'J written as 

Ip"' [ s "] ip 

- Pgt'^o) R2("i^l) (B'2) 
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Fig, 3.6.9 Simulation Errors in Rigid Body Case 
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Fig. 3.7.1 Simulation Errors in Rigid Body Case 


-6 


RUN NO. 2684 


6St 

EULER PARAMETERS ERRORS 



RUN NO 


2 684 













f9X 

fiNGULflR VELOCITY ERRORS 110"'^ RflO/SECl 






991 

RNGULflR VELOCITY ERRORS (1 0""^ RflD/SECl 




EULER PARAMETERS ERRORS 

MTS ASYMPTOTIC SOLUTION 


DT = 500 sec jg q _ 

GEOMAGNETIC TORQUE 0 j - 




6 


PNGULflR VhLUCITT ERf=?QRS 


t 


S f 


C_1 

LU 

in 


Q 

CC 

OC 


o 


4 -- 

3 

2 


TORQUE FREE SOLUTION 
GEOMAGNETIC TORQUE 


WX = 

WY ■= -i 1 K 

WZ “ Q © o 


CTi 

CO 


in 

ac 

Q 

QC 

az 

LU 

>- 


o 

Q 

_J 

LU 

> 

CC 

CC 

_1 

o 


0 


-2 + 


-3 4 - 



-4 

-5 


Fig. 3.7,11 Simulation Errors in Rigid Body Case 
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Fig. 4.4.8 Simulation Errors in Dual Spin Case 


RUN N0.= 2881 



EULER PARAMETERS ERRORS 


EULER PARAMETERS ERRORS 


Dual -Spia Satellite 

Torque Free Solution - Reference 

Gravity Gradient Torque 


- 

- 



TIME (SEC) 


Fig, 4.4,10 Simulation Errors in Dual Spin Case 


RUM M0.= 


081 

RNGULflR VELOCITY ERRORS (lo" IRRD^SECl 


6 


TlMGULqiR VEECCITY ERRORS 


5 

4 

3 

2 


Dual Spin Satellite 

MTS Asymptotic Solution " Reference 
D T = 10 sec 

Geomagnetic Torque 






EULEf=i PARAMETERS ERRORS 


Dual Spin Satellite 

MTS Asymptotic Solution - Reference 
D T= 10 sec 

Geomagnetic Torque 




•Fig, 4,4.12 Simulation Errors in Dual Spin Case 
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Table 3.6.1 Rigid Body Satellite Disturbed By Gravity 
Gradient Torque 


Runs Description 



Direct Integra- 
tion 

^ T = 10 sec 
(reference case) 


Asymptotic 

Approach 

AT = 10 sec 


Direct Integra- 
tion 

A T = 50 sec 


Asymptotic 

Approach 

A T = 100 sec 


Asymptotic 

Approach 

A T = 200 sec 


Asymptotic 

Approach 

A T = 500 sec 


Torque-free Sol. 
Subtract Ref. 
Case (Effect of 
G.G.T.) 


* IBM 360/75 run time. 

** Does not include initialization time - for nominal solu- 
tion and Fourier transformations - in this case is about 
40 sec. 


201 

































Table 3.7,1 Rigid Body Satellite Disturbed By 
Geomagnetic Torgue 


Runs 

( 

Description 


■m 

trrors I 


1 

Direct Integra- 
tion 

. 4 T = 10 sec 
(reference case) 

112 

0 

0 


2 

Asymptotic 

Approach 

AT = 10 sec 

294 

'7 

0.2 X 10 

.002 

BB 

BBS 

3 

Direct Integra- 
tion 

At = 25 sec 

45 

3.0 X lo'"^’ 

.005 

n 

m 

4 

Asymptotic 

Approach 

A T = 100 sec 

31 

0.2 X 10 

.002 


5 

Asymptotic 

Approach 

A T = 200 sec 

16 

-1 

.25 X 10 

.003 


6 

Asymptotic 

Approach 

A T = 500 sec 

6.5 

.30 X 10 

.003- 

B 

McwmwnI 

7 

Torcrue-free Sol, 
Subtract Ref . 
Cass (Effect of 
G.M.T.) 


100 X 10 

.003 

m 
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Table 4.4.1 Dual Spin Satellite Disturbed by Gravity 
Gradient Torque 


Runs Description 


Direct Integra- 
tion 109 

A T = 10 sec 

(reference case) 


Asymptotic 

Approach 

A T = 10 sec 


Direct Integra- | 

3 tion i . 

i 22 10 X 10 

T = 50 sec ! 


Max. Errors 




Torque-free Sol, 
Subtract Ref. 
Case (Effect of 
G.G.T.) 


4.4.2 


4.4.3 

4.4.4 


4.4.5 


0.3 X 10 .006 4.4.6 



,5 X 10 


.5 X 10 


4.4.7 

4.4.8 


FtCr. 

4.4.9 


4.4.10 
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Table 4.4.2 Dual Spin Satellite Disturbed By 
Geomagnetic Torque 


Runs 

Description 

Computer 

Time 

(sec) 

Max. Errors 

Figures 



1 

Direct Integra- 
tion 

Zk T = 10 sec 
.(reference case 

102 

> 

1 

0 

0 


2 

Asymptotic 

Approach 

A T = 10 sec 

290 

-S' 

0.5 X 10 

.0010 

Fig, 

4.4.11 

4.4.12 

3 

Direct Integra- 
tion 

A T = 25 sec 

42 

LO X 10 

.0150 

Pig. 

4.4.13 

4.4.14 

4 

Asymptotic 

Approach 

A T = 100 sec 

30 

-S 

0.5 X 10 

.0015 

Fig. 

4.4.15 

4.4.16 

5 

Asymptotic 

Approach 

A T = 500 sec 

6.4 

-i- 

0.6 X 10 

.0017 

Fig. 

4.4.17 

4.4.18 


204 






























Table 5.5 GMS Approximation*and Direct Solution 


Case 

Description 

Computer 
time (sec) 

Fig. 


GMS Approximation 
/it=1000 sec 

2.9 

5.5.3 

5.5.4 

2 

GMS Approximation 
At=2000 sec 

1.6 

5.5.5 

5.5.6 


Direct integration 
A t=10 sec 

57.9 

5.4.3 

5.4.4 

4 

Direct integration 
At=30 sec 

19,2 

5.5.7 

5.5.8 


* General multiple scales with nonlinear clock. 
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